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effective  wind  vector  and  shear  parameter. 
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UNCLASSIFIED 


SECURITY  CLASSIFICATION  OF  THIS  P  AGE'K/II-"  ».>'*  t " ri-rrrf) 


PREFACE 


The  author  gratefully  acknowledges  the  support  and  cooperation  of 
Dr.  David  L.  Auton  of  the  Defense  Nuclear  Agency,  and  Mr.  Ralph  B.  Mason 
of  the  Command  and  Control  Technical  Center  (CCTC)  who  made  exhaustive 
runs  of  the  code  on  the  CCTC  computer,  analyzed  the  results,  uncovered 
numerous  problems,  and  thus  substantially  assisted  in  the  development 
of  the  model . 


TABLE  OF  CONTENTS 


Section 


1. 

2. 


3. 


4. 


5. 


PREFACE - - - 

INTRODUCTION  AND  BACKGROUND  -  -  -  - 

DATA  BASE  - 

2.1  ACTIVITY  DEPOSITION  RATE  -  -  -  - 

2.2  FALLOUT  ONSET  TIME - - - 

2.3  FIREBALL  AND  STABILIZED  CLOUD  DATA  - 

2.4  NOMINAL  PARTICLE  -  - 

2.5  OBSERVED  FALLOUT  DATA  - 

MATHEMATICAL  DESCRIPTION  - 

3.1  ACTIVITY  DEPOSITION  RATE  FUNCTION  - 

3.2  DEPOSITION  AS  A  FUNCTION  OF  DISTANCE  FROM  GROUND 

ZERO  -  - - - 

3.3  FARF1ELD  CORRECTION - -  - - 

3.4  UPWIND  CORRECTION  -  -  - 

3.5  FALLOUT  TIME  OF  ARRIVAL  - 

3.6  HORIZONTAL  SPREAD  OF  THE  NUCLEAR  CLOUD  - 

3.7  TURBULENT  DISPERSION  OF  FALLOUT  -  -  -  - 

3.8  WIND  SHEAR  DISPERSION  AND  CROSSWIND  SPREAD  OF  THE 

FALLOUT  PATTERN  - 

3.9  GAMMA  RAY  EXPOSURE  RATE  AND  MAXIMUM  EFFECTIVE 

BIOLOGICAL  DOSE - - - 

USE  OF  WIND  DATA - -  -  -  - 

4.1  GENERAL  CONSIDERATIONS  - 

4.2  FEFECTIVE  FALLOUT  WIND  - 

4.3  SHEAR  PARAMETER - - - 

VALIDATION - - - - 

5.1  DISCUSSION  OF  RESULTS - - - - - 

5.2  DISCUSSION  OF  THE  TEST  SHOT  DATA  AND  PREDICTIONS 


1 

7 

9 

9 

9 

19 

20 
22 
23 
23 

25 

30 

34 

35 
35 

40 

41 

44 

47 

47 

47 

50 

52 

52 


2 


TABLE  OF  CONTENTS  (Continued) 


Section  Page 

5.3  OBSERVED  AND  PREDICTED  FALLOUT  PATTERNS  -  57 

6.  COMPUTER  CODE - -  81 

6.1  GENERAL  DISCUSSION - 81 

6.2  FALLOUT  MAPS . . . 83 

6.3  INPUT  OF  WIND  PROFILE  DATA  -----  - -  85 

6.4  STORAGE  AND  COMPUTATION  TIME  REQUIREMENTS  -  90 

6.5  DESCRIPTION  OF  CARD  INPUTS  FOR  THE  DNAF-1  CODE  -  -  91 

6.6  EXAMPLE  PROBLEM  AND  PRINTOUT  DESCRIPTION  -  93 

REFERENCES - 97 

APPENDIX  A  GLOSSARY  OF  SYMBOLS  AND  FORTRAN  MNEMONICS  -  .  -  99 

APPENDIX  B  GROUND  ROUGHNESS  AND  INSTRUMENT  RESPONSE  CORRECTION 

FACTORS - 103 

APPENDIX  C  FALLOUT  PATTERN  COMPARISON  BY  THE  FIGURE-OF-MERIT  METHOD  105 

APPENDIX  D  FORTRAN  CODE  FOR  THE  DNAF-1  FALLOUT  MODEL  - -  107 


LIST  OF  ILLUSTRATIONS 


j 


•  , 


« 


* 


4 


i 


;<% 


Figure  Page 


1  Computer  plots  of  activity  fraction  deposition  rate 

vs.  time  as  computed  by  DELFIC.  -----------  10 

2  Activity  fraction  deposition  rate  function,  g(t), 

(eq.  (6)),  without  farfield  correction,  vs.  time 
for  W  =  1  KT.  A  sampling  of  DELFIC  results  are  in¬ 
cluded  for  comparison.  ---------------  24 

3  Distribution  of  activity  of  depositing  fallout  in 

the  hotline  axis  direction.  -------------  27 


4  Comparison  of  Gaussian  function  with  the  function 
used  in  DNAF-1  to  approximate  the  spatial  distribu¬ 
tion  of  cloud  activity  along  the  hotline  axis.  -  -  - 

5  Activity  deposition  rate,  g(t)f,  including  farfield 
correction,  vs.  time.  A  sampling  of  DELFIC  results 


are  included  for  comparison.  ------------  32 

6  Basis  of  the  timc-of-arrival  calculation.  ------  36 

7  Fallout  time  of  arrival  vs.  distance  from  ground  zero 

for  several  yields  as  computed  by  DNAF-1  and  WSLG-10.  37 

8  Crosswind-integrated  activity  fraction,  D(X)*-  and 
D(X)U,  vs.  distance  from  ground  zero  along  trie 

hotline.  ----------------------  42 

9  Wind  data  card  input  for  the  example  data  listed  in 

Table  7. - - - - -  89 

10  Test  problem  card  input.  --------------  94 


•  ; 


j 


•A 


i  g 


i 


i 


«  4 


■  1 
\ 

•vi 

a 

..  i 
i 

i 


LIST  OF  TABLES 


Tab!  e  Page 

1  SETTLING  SPEEDS  FOR  THE  NOMINAL  PARTICLE,  & 

-  229  um - -  - -  21 

2  TEST  SHOT  DATA  -  52 

3  COMPARISON  OF  OBSERVED  AND  PREDICTED  FALLOUT  JATTERN  54 

STATISTICS  - - - 

4  OVERALL  MEAN  ABSOLUTE  PERCENT  ERRORS  -  55 

5  EFFECTIVE  FALLOUT  WINDS  AND  SHEAR  PARAMETERS  COM¬ 
PUTED  FROM  H  HOUR  WIND  PROFILES  FOR  USE  BY  DNAF-1  -  57 

6  DESCRIPTION  OF  DNAF-1  CODE  SUBROUTINES  AND  FUNCTIONS  82 

7  EXAMPLE  WIND  DATA  LISTING -  88 


1  .  INTRODUCTION  AND  BACKGROUND 


To  quickly  and  efficiently  estimate  fallout  radioactivity  from  large 
numbers  of  nuclear  surface  explosions,  for  example,  for  military  damage 
assessment  studies,  a  simple,  very  fast  fallout  prediction  code  is  needed. 
While  codes  based  on  numerical  models  ’  ’  provide  flexibility  of  usage 
and  relatively  high  prediction  accuracy,  they  are  cumbersome,  use  too  much 
computer  storage,  require  more  input  data  than  desired,  and  use  too  much 
computer  time  per  prediction.  A  model  which  uses  analytical  equations 
rather  than  a  numerical  approach  is  appropriate  for  this  purpose. 

The  model  that  has  best  satisfied  these  requirements  in  the  past,  the 
WSEG-10  model ,  has  been  used  for  more  than  twenty  years  for  damage  assess¬ 
ment  studies4’5.  WSEG- 10  has  recently  been  analyzed6  and  its  prediction 

7 

capabilities  compared  with  those  of  several  other  models  .  It  was  found 
that,  while  in  several  respects  WSEG- 10  is  satisfactory  in  terms  of  its 
mathematical  structure,  its  data  base  is  obsolete,  and  this  deficiency 

7 

alone  was  seen  to  substantially  compromise  its  prediction  capability  . 

To  upgrade  prediction  capability  the  easiest  course  would  be  to  upgrade 
the  WSEG-10  data  base,  but  retain  its  mathematical  structure.  However 
owing  to  several  deficiencies  of  the  mode1  itself,  this  course  has  turned 
out  to  be  undesirable.  The  most  important  of  these  deficiencies  are  as 
follows.  WSEG-10  mathematics  are  based  ori  a  curve  fit.  of  an  exponential 
function  to  radioactivity  deposition  rate  data  that  were  calculated  by  an 
early  fallout  model.  The  particular  model  used  was  developed  to  predict 
fallout  iron  surface  bursts  of  large  yield  (i.e.,  megaton  range)  nuclear 
weapons.  Tor  these  large  yield  cases  the  exponential  function  fits  the 
deposition  rate  data  reasonabl y  wel 1  at  intermediate  and  late  times.  Un¬ 
fortunately,  this  is  not  the  case  for  low  yield  explosions.  Also,  in  no 
case  does  the  exponential  function  go  to  zero  at  detonation  time  as  a 
physically  realistic  function  should  do.  Thus,  to  fit  the  low  yield  deposi¬ 
tion  rate  data  as  well  as  the  high  yield  data,  and  to  force  the  fitted 
function  to  go  to  zero  at  zero  time,  a  new  mathematical  base  is  needed. 


Another  serious  problem  with  WSEG-10  is  that  it  assumes  that  a  Gaussian 
function  describes  the  vertical  distribution  of  activity  in  the  nuclear 
cloud.  For  shots  with  yields  less  than  about  50  KT ,  this  is  a  very  poor 
assumption  since  much  of  the  activity  which  will  fall  out  locally  is  in  the 
cloud  stem,  whereas  the  WSFG-10  Gaussian  peaks  near  the  cloud  cap  center 
height.  In  effect  this  procedure  ignores  the  stem,  even  for  low  yield 
cases,  and  this  also  seriously  compromises  prediction  capability. 

To  substantially  improve  prediction  capability,  and  to  extend  the  range 
of  applicability  to  1  ov/er  yields,  an  entirely  new  model  is  required.  The 
DNAF-1  model  has  been  developed  to  fill  this  need. 

DNAF-1  also  is  based  on  a  curve  fit  to  calculated  activity  deposition 
rate  data,  though  a  new  set  of  data  is  used,  and  a  more  appropriate  function 
is  used  to  fit  the  data.  An  entirely  new  approach  is  used  to  account  for 
vertical  structure  of  the  stabilized  cloud,  which  does  account  for  fallout 
from  the  cloud  stem.  Indeed,  revised  data  and  updated  modeling  concepts 
are  used  throughout  the  development  of  the  new  model.  Also,  the  yield 
range  has  been  extended  downward  from  1  KT,  the  lower  limit  for  WSEG-10, 
to  10“3  KT.  The  upper  yield  limit  is  10J  KT. 

The  model  computes  gamma  radiation  from  dry  ground-deposited  particu¬ 
late  fallout  from  the  nuclear  cloud  cap  and  stem.  This  radiation  issues 
from  weapon  debris  fission  products  and  from  induced  activity  in  the 
fallout.  Contributions  from  throwout  and  induced  activity  in  the  crater 
or  elsewhere  are  not  included.  Alpha  and  beta  radiation  are  not  treated. 

In  the  next  section,  section  2,  we  discuss  the  more  critical  aspects 
of  the  data  base  used  for  the  new  model.  In  section  3  we  detail  the  mathe¬ 
matical  structure  of  the  model,  and  in  section  4  we  discuss  how  wind  data 
are  used.  Section  5  contains  results  of  a  validation  study:  predicted 
fallout  patterns  are  compared  with  observed  patterns  for  five  test  shots, 
and  predictions  by  the  WSEG-10  and  ()ELrIC  models  are  included  as  well  as 
those  by  the  DNAF-1  model.  Finally,  the  computer  code  is  described  in 
section  6. 

Throughout  the  presentation,  we  make  comparisons  between  DNAF-1  and 
WSEG-10  modeling  approaches  arid  results  wherever  such  is  appropriate. 


2.  DATA  BASE 


2.1  ACTIVITY  DEPOSITION  RATE 

As  mentioned  above,  the  DNAF-1  model  is  based  on  a  mathematical  function 
that  is  fitted  to  activity  deposition  rate  data.  These  data  represent  rate 
of  deposition  of  total  fallout  activity  as  a  function  of  time  and  yield.  Of 
course,  it  would  be  best  to  use  observed  data,  but  there  are  none,  so  instead 
we  have  used  results  generated  by  a  special  version  of  the  DELFIC  code1’2. 

Nuclear  cloud  rise  was  calculated  through  the  1976  U.S.  Standard  Atmos¬ 
phere8  for  every  decade  of  yield  from  I0’3  KT  to  IQ5  KT.  Explosion  was 
taken  to  be  at  sea  level.  Fission  yie'M  was  taken  to  equal  energy  yield. 
Twenty  thousand  fallout  parcels  (100  particle  size  classes  and  200  cloud 
subdivisions)  were  followed  for  each  calculation.  Cumulative  activity 
tabulated  as  a  function  of  time  was  differentiated  to  give  deposition 

g 

■rate  by  the  method  of  cubic  splines  .  Results  are  given  in  Figure  1. 

The  uneven  appearance  of  the  curves  at  early  times  was  caused  by  poor 
deposition  statistics  for  the  sparse,  though  highly  radioactive,  early 
fallout.  Some  numerical  experimentation  showed  that  increased  resolution 
of  the  calculations  produced  curves  of  essentially  the  same  shape;  there¬ 
fore,  to  save  computer  time  the  curves  shown  were  accepted  as  an  adequate 
compromi se. 

The  fitting  of  a  mathematical  function  to  these  data  is  described  in 
section  3.1. 

2.2  FALLOUT  ONSET  TIME 

Fallout  onset  time,  tQ,  is  the  time  of  first  impact  of  fallout  on  the 
ground.  It  is  needed  to  calculate  fallout  time  of  arrival  as  a  function  of 
distance  from  ground  zero,  which  in  turn,  is  required  For  calculation  of 
several  critical  parameters.  Onset  times  for  each  decade  of  yield  were 
tabulated  from  the  DELFIC  deposition  rate  data,  and  n  (t(i)  was  fitted  by 
least  squares  to  a  polynomial  in  tn(W).  The  result  is 
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ricjure  1.  Computer  plots  of  activity  fraction  deposition  rate  vs.  time  as  computed  by  DELFIC.  '' 
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icjure  1.  Computer  plots  of  activity  fraction  deposition  rate  vs.  Li  mo  as  computed  by  DllflC. 
(continued ) 


» 


11 


f 


■SIElO- 100 -OQO  KT 


1.12  i.S3  I. *4  2.36  2-77  3. 1 8  J.59  4.00  4.4|  4.R> 


-6.29 


t  ■■  GXp 


1.52/667  '■  0. 40S94(>6-  nW  '  0. 02064  322(  11W ) 


;  W  10  K1 


(1) 


t  -  1147.54 


;  W  10 


wlmre  t  is  in  units  of  seconds. 

2.3  TIREBALL  AMD  STAI5 II. l/TD  ('I  'U|:  DATA 

The  model  a  t  •'  o  requites  fireball  (i.o,,  early  cloud)  radius,  as  well  us 
radius  and  lop  and  ha  ,e  lieiqht  .  of  !  lie  stabilized  cloud.  The  latter  are 
tukt  n  to  be  above  mean  sea  level. 

The  early  cloud  radius.  K .  (m) ,  is  computed  from  a  function  used  to 
initialize  the  S1W  1C  ol oiul  rise', 

1^  --  108Wn’  '  W  in  KT.  (?) 

Stabilized  cloud  ran  boa;  lieiqht,  ,  and  top  lieiqht,  /y,  are  computed 
from  functions  fitted  to  observed  data  by  Wilsoy  and  Crisco10, 


zu 

=  aW 

(3) 

zi 

,.d 
-•  cW 

(4) 

whore  W  is  energy  yield  in  kilotons. 

the  heights  are  in  meters,  and 

a  = 

2228 , 

b  -  0.3463; 
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— I 
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2661  . 
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Stabi  1  i  >: 

:ed  cloud  radius,  ltc. , 

i  s 

computed  from  a  least  squares  poly- 

nomial  curve 

fit  to 

Dill  1-  1C  result: 
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?  1 

R  =  exp  6.7553  +  0.32055'^nW  +  0. 01 1 37478  UnW)  '  I  (5) 

where  is  in  meters  and  W  is  yield  in  kilotons. 

2.4  NOMINAL  PARTICLE 

The  model  requires  that  a  single,  "effective  fallout  wind"  vector, 
v,  be  calculated  from  a  given  vertical  profile  of  wind  data.  This  effective 
fallout  wind  vector  is  an  average  of  the  given  data,  where  in  the  averaging 
each  vector  is  weighted  according  to  the  time  required  for  a  "nominal  particle" 
to  settle  through  the  wind  space  stratum  represented  by  that  vector. 

The  nominal  particle  diameter  was  calculated  from  DELFIC  results  to  be 
the  activity-weighted  average  diameter  over  all  particle  sizes.  The  standard 
DELFIC  lognormal  particle  size  distribution  was  used, 

i*50  =  0. 407  mu 

s  =4.0 

where  6<>0  is  the  median  diameter  and  s  the  geometric  standard  deviation  of 
the  distribution  of  particle  number  with  respect  to  particle  diameter.  (See 
Appendix  A  of  reference  1.) 

Calculations  were  performed  for  the  full  range  of  yields,  at  every 
dn  ade  of  W  starting  with  10"'5  KT ,  using  a  special  version  of  DELFIC. 

The  low  yield  results  were  averaged  separately  from  the  high  y;eld 
results,  to  give 

&  =  217.05  Mill ;  W  <  10C  KT 

nom  ’  — 

&  --  240.  95  uni;  W  >  100 

nom 

The  average  of  these  diameters  is  229  i.m,  which  is  the  value  used  in  the 
model  . 
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Using  6  =  229  mil,  a  table 

s  nom  ’ 

as  a  function  of  altitude  m  the 

(Table 

1976  U. 

1)  of  particle  settling 

g 

S.  Standard  Atmosphere  , 

speeds 

was  com- 

puted  via  the 

equations  given  in 

section 

2.2.3  of  the  DELF1C  documentation 

This  table  is 

used  by  subroutine 

EFWIND 

to  calculate  effective  fallout 

wi  nd . 

TABLE  1 

SETTLING  SPEEDS  FOR  THE 

NOMINAL 

PARTICLE,  6  =  229  pin 

A1 titude 
(km  above  MSL 

Sett 1 ing 
Speed 

l  (m  s-1) 

A1 ti tude 
( km  above  MSL) 

Settl ing 
Speed 
(in  s-1 ) 

0 

1.6538 

20 

3.5222 

1 

1.7124 

21 

3.6322 

2 

1.7744 

22 

3.8122 

3 

1.8401 

23 

3.9243 

4 

1.9097 

24 

3.9994 

5 

1.9836 

25 

4.0806 

6 

2.0621 

26 

4.1746 

-» 

/ 

2.1458 

28 

4.3565 

8 

2.2350 

30 

4.5328 

9 

2.33U3 

32 

4.7083 

10 

2.4324 

34 

4.8517 

11 

2.5419 

35 

5.0023 

12 

2.6446 

38 

5.1767 

13 

2.7489 

40 

5.3910 

14 

2.8551 

42 

5.6632 

15 

2.9630 

44 

6.0155 

16 

3.0723 

46 

6.4727 

17 

3. 1831 

48 

7.0778 

18 

3.2951 

50 

7.8819 

19 

3.0482 

OBSLIiVID  FALLOUT  DATA 
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After  developmrn I  of  the'  model  ,  its  capah  i  1  i  I  y  must  !>o  cv.i  1  ua  I  !  hy 
comparison  anainst.  observed  data.  lor  this  purpose,  we  have  used  (he  host 
who  1 1'- po  f  t.c'tn  fallout  data  available,  Comparison  results  are  not  qroaliy 
different,  fiom  those  obtained  hy  use  ol  presumably  much  more  capable 
model  s  ’  ’  .  Details  are  p  i  von  below  in  sect  ion  !>. 
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i.  MAI  H[' MAT  i  FAI  ill  S'  is  f  I ' !  I  ON 

i  i  act  i  v i  r y  v-  pus  i  r  ion  rah  i  i;nc :  1 

[tloa  1  a-ec  i  f  i  ca  t.  i  ons  lot  ( l,r  aoijviiy  d'po-.  i  i  ion  into  funcfion,  i  j  ( t. )  , 
<ii!'  .is  fnl  1  ovr,  : 

I.  If  .  lieu  Id  bp  p  ,  imp  Ip  single  lunrt  ion  ,  ppi|ii  i  r'i  ric|  .is  fow  yield 
*1 1  •  |  n'ndt't  1 1  parameters  as  pram  icahl  which  fits  the  deposition 
rail-  da !  a  (rip.  I)  over  the  mUrr  niiip  ilpcados  of  yield  ramie. 

I  lie  product  of  'i(l)  with  a  downwind  di  ,  pens  ion  factor,  F(X,t. ), 
should  be  analytically  i  nteur  ib  I  e  over  tilin'  to  yield  activity 
a.  a  funcfion  of  downwind  distance.  (See  the  next  section.) 

a.  bo i h  downwind  and  crosswind  dispersion  functions,  F  ( X  , t )  and 
(iff),  should  have  (V.usSian  forms  as  required  by  theory  and  con¬ 
firmed  by  experiment  for  diffusion  by  homogeneous  turbulence11, 
in  the  absence  of  firm  data  to  support,  alternative  functions. 

Thus,  specification  2  is  extended  t.o  require  that  q(f  )  be  such  as 
to  allow  analytical  integration  of  the  product  g(t)F(X,t)  where 
I  ( X  ,  t, )  is  spec  i  f  jcal  lv  a  Gaussian  function. 

4.  u(t)  should  go  to  zero  at  f. ho  limits  of  zero  and  infinite  time. 

!>.  Not  mal  ization  should  be  exact  such  that  activity  is  conserved, 
these  specifications  are  so  restrictive  as  to  require  some  compromise.  Fven 
■■o,  only  i'"  jti  i  reutonf.  s  3  and  b  are  not.  met  by  the  Function  selected. 

ttu;  function  which  comes  closest  to  satisfying  the  specifications 
was  found  to  be 

4.  (ft) 

lien?  g(t)  is  t  ate  of  d>  position  of  activity  fraction  (s’-1),  t  is  time  (s), 
and  <  (dimensionless)  arid  I  (s)  (iro  yield  dependent  constants. 

Injure  3  shows  a  representative  plot  of  eq.  (6),  for  W  1  KT ,  and 
compares  this  with  a  samplinq  of  the  1)11  I  IC  results.  FxcepL  for  a  tendency 
to  overprt'd  i  ct  the  m,,ximum  deposition  rate  (also  see  fiq.  b )  For  some  cases, 
the  function  i  >,  capable  of  reproduc  i  ng  the  1)11  Ilf  results  reasonably  well 
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DELFIC  results  are  included  for  comparison. 


for  early  and  intermediate  times.  At  late  time  the  eq.  (6)  function  becomes 
linear  (on  the  log-log  plot)  whereas  the  deposition  rate  data  drop  off  rapidly 
This  discrepancy  requires  application  of  a  "far  field"  correction  factor 
which  is  discussed  in  section  3.3. 

The  maximum  of  the  g(t)  function  occurs  at  time  t  which  is  qiven  by 
3  max  3  J 

tax  =  W4ir  (O 

The  significance  of  this  relation  is  that  only  two  of  these  three  parameters 
need  to  be  specified  as  functions  of  yield.  Approximate  best  fit  values  of 
a  and  t  were  determined  somewhat  subjectively  by  trial  and  error  Tor  each 
decade  of  yield.  These  were  fitted  to  simple  functions  of  yield  as  follows: 


a  =  1 

.06  ; 

10 

<_  w  < 

10"'  KT 

a  =  1, 

.0875  +  .0119431  £nW; 

10"1 

<  w  < 

103 

(8) 

n  1. 

■  17 

10 3 

<  w 

and 


t  - 

max 

0.41 SSG 

30W  ; 

io- 1 

•:  w 

< 

1  KT 

t  = 

max 

0*0440/ 

30W 

1 

■:  U 

< 

io3 

t 

max 

0.  1 G  ?  7  0 

893.616W 

io3 

w 

io“ 

(9) 

^max 

0  •  0  b  1 1  s 

2497. 18W 

ID3 

W 

10S 

where  t  is  in  units  of  seconds, 
ma  x 

3.2  DEPOSITION  AS  A  FUNCTION  OF  DISTANCE'  FROM  GROUND  ZERO 

The  action  of  effective  fallout  wind  v  is  to  transport  the  stabilized 
cloud  downwind  such  that  at  time  t  it  is  centered  over  downwind  distance 
point  X  ---  vt.  We  assume  that  fallout  deposited  from  the  cloud  at  this  time 


?b 


I 
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...  4 


>1 


has  a  continuous  distribution  along  (Ik  X  (i.e.  hotline)  axis  with  the 
distribution  centered  and  peaked  at  K  t  .  Hie  situation  is  illustrated 
in  figure  3.  To  determine  the  total,  crosswind  integrated  fallout 
deposited  )L  any  distance  X  from  ground  'em,  j.  ho  coniniuil  ions  of  .ill 
fallout  deposits  must  be  summed  over  ail  I  imr,  and  in  this  case-  this  r, 
accomplished  by  analytical  integration  ‘ ■ om  t  -  0  to  t  -  ■  rf  the  proa  act 
of  q(t)  with  the  downwind  d  istr  ihiiti  ' n  fun*  'Lion. 

As  explained  in  the  previous  soi  l  ion,  a  most  desirable  selection 
for  the  downwind  distribution  function  for  depositing  fallout  would  be  a 
Gaussian  function.  It  turned  out,  however,  that,  we  could  not  find  a 
satisfactory  function  for  g(t)  that  is  analytically  integral)  1  e  in  combina¬ 
tion  with  a  Gaussian  distribution  in  X-vt.  Thus,  we  were  forced  to  use 
the  distribution  function 


r  ( x ,  t ) 


no  1  + 


(10) 


where  a  is  analogous  to  the  Gaussian  standard  deviation.  Figure  4  shows 
a  comparison  of  the  eg.  10)  function  with  the  Gaussian  function. 

Thus,  total,  crosswind  integrated  activity  fraction,  D(X),  deposited 
at  distance  X  along  the  windward,  or  hotline,  axis  is 


A  result  is  obtained  in  closed  form  provided  we  take  «=1  inside  the  integral. 
This  result  is 
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NAF-1  to  approxi./ia 


.M 


v2xl  (3-«)  | 


-  ttXX !  (3Xp;  +  X?)  X2  +  4X2X] 


X3  +  4X2Xj)  X2  +  (xfx$  -  4 X 2 X j  )  £n  A 


~  _(X2+3X?  )  X2X3  +  4X2X';_  11 


+  2  tan 


(12a) 


where 


Xj  =  vT 

2  2  2 

X2  =  X  4  O 

2  ?  2 

X3  =  Xz  -  Xj 

x5  =  (xl3  4  8XZX? )x5  +  1 6 X 4 X  j  , 
and  D(X)  has  units  in”1. 

This  is  the  basic  equation  for  the  DNAF-1  model.  In  the  computer  code 
(function  DNAF1),  an  alternative  form  of  eq.  (12a)  is  also  used,  which  is 
required  because  of  limitations  of  some  computers 

D( X )  ~  -2iLiJP  M^l/ ij  |  _,rXXl  (3X2  4  X?)/xJ  4  4X?Xi/x'; 

*  Xs(3-«)  L  J 

1  /  o  \ 


where 


?  /  72 

X  (1  4  4X  Xi 


2/X4  )  4  (xz  -  4X2Xj /X3) 


(12b) 


■  9 V2  \ /  V2  +  /IV2V!!  /v4  „  +  9  |  >n"  1  /X  \1  '  .  I  V2  I 

^  A  2  1  3  A  ]  J  j  /\  j  +  1  /\  i/  /x  q  7T  »  t.  Lull  ^  ^  yj  ^  9  V  'v  ^  I 


4  4  ?.  ?  4  4  4 

Xb  =  X;!  4  8X  Xi  4  16X  X,/X:, 
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Thus  eq.  (l?a)  is  used  for  small  values  of  |X(,|  and  eq.  (12b)  for  larqe 
values. 

Negative  values  of  X  (i„o.,  upwind  distances)  are  accomodated  as  well 
as  positive  values  in  eqs.  (12a)  and  (12b). 

For  a  point  cloud,  defined  by  o  -  0,  eq.  (12)  reduces  to 


— -.-X-  -  ;  x  -  o 

(X'fX',)‘ 


(13) 


which,  on  substitution  of  X  -  vt,  becomes 


vl)(vt)  , 

cl  —  U 


4s_i  n  [  iij  3 -.»)/?] 
"  (3->0 


(t  iT  )' 


(1/1) 


The  riqht  hand  side  of  eq.  (14)  is  equivalent  to  eq.  (6)  with  (T/t)u  replaced 
by  T/t.*  Thus,  provided  that  ■*  is  not  much  different  from  unity,  eqs.  (13) 
or  (14)  may  be  used  instead  of  eq.  (6)  to  represent  activity  deposition  rate 
as  a  function  of  either  distance  from  ground  zero  or  time  for  a  point  cloud. 
These  equations  are  used  below  for  the  development  of  the  farficld  correction. 

3.3  FARFIELD  CORRECTION 

As  already  noted,  the  eq.  (6)  function  for  deposition  rate  fits  the 
DELFIC  results  adequately  at  early  and  intermediate  times,  bul  at  late  times 
(which  correspond  to  farfield  deposition)  the  DELFIC  results  drop  off  much 
more  rapidly.  Moreover,  this  discrepency  becomes  more  acute  as  yield  increases 
Indeed,  as  inspection  of  eq.  (14)  shows,  g  ( t )  *  1/t;  t  ‘  ,  whereas  the 

DELFIC  result.-  for  largo  time  are  proportional  to  oxn  [-(t./T)n~|,  n  =  l  or  ? 
and  i  constant. 

A  correction  to  the  functions  for  y(t)  and  D(X)  which  have  general 
applicability,  but  are  effective  only  for  late  times  and  correspondingly 
large  X  are  as  follows.  As  shown  at  the  end  of  the  last  section, 


*Recall  that  >t  was  taken  to  bo  unity  inside,  the  integral  iri  eq.  (11). 


n ( X )  .q,  (i.e.,  downwind  deposition  for  a  point  cloud)  is  oqui va  1  ent 

to  y(t)/v.  It  ill  so  turns  out  that  D(vL)  ,.  q(t)/v  for  la 

o  r  vJ 


i  ■  y  X  -  v  t . 

The  correction  factor  is  derived  by  moms  of  an  exponential  interpola¬ 
tion,  1  -exp  [  -X/ (avt(. }]  ,  between  Lhe  functions  n 

where  a  is  constant  and  tr  ,  k  and  x  are  function''  of  yield.  The  farfield  ■ 

- 1 


■  t.  1 

Chus , 

[rue  t 

,  which 

t.  are 

rela I er 

a  1  ini 

lerpola- 

1  k 

r  /X  \ 2 ] 

vxp 

L"U  J 

orrected  ll(X),  1J ( X ) ^  (m"  ),  is 

n  (  vl)(X )  / k  )  -  (*-.) 


D(X)f  s  |j(x)exp 


1 -exp  |  -X/ |avtc J 


X  >  0 


(1! 


computed  by  eg.  (12a) 

or  (12b) 

1.443 

-n  ')[  -n,  r; 

9.867X10  W  • 

(s'1) 

O 

O 

<o 

(s) 

14 667W0’ 26208 

( s ) ;  W 

exp  10.124706  t  0. 1861  768sinW 
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-  O.U08660444  (•••■» W)' 


(s)  ;  W  .  98.787  KT 


Figure  5  shows  he  function  <j(t),.  (actually  vD(vt)J.  see  eci. 

•  ‘  1  » f  T  “  ( J  1 

14))  for  every  other  decade  of  yield  along  with  a  sampling  of  the  I1FI.FIC 
results  and  the  corresponding  quantity  for  the  WS1G-10  model. 

The  WSFG-10  model  uses  the  deposition  rate  function 


9  ( t ) 


1 


WSFG-10  "  TTT1 +T7n~) 


exp 


[-  (t/T')n 


(17 


where  I  is  the  gamma  function,  T‘  is  a  yield  dependent  constant  and  n 
has  a  value  between  1  and  2.  In  calculate  g  the  results  shown  in  figure 
5,  we  have  used  n  -  1.6.  Notice  in  Figure  b  that,  the  WSFG-10  function 


matches  the  DELFIC  deposition  rate  data  only  in  the  farfield  region,  and 
that  instead  of  going  to  zero  at  t  =  0,  it  actually  peaks  there.  Also 
note  that  the  WSEG-10  curve  fit  probably  is  adequate  for  very  high  yields, 
but  is  very  poor  for  low  yields. 


3.4  UPWIND  CORRECTION 


In  comparing  calculated  fallout  patterns  with  observed  patterns  for 
low  yield  test  shots,  it  was  found  that  the  calculated  activities  de¬ 
creased  too  slowly  with  upwind  distance  from  ground  zero.  Analysis  of  the 
discrepancies  indicated  a  correction  factor  that  is  independent  of  both 
yield  and  wind  speed  for  yields  less  than  10KT.* 

For  high  yield  shots  there  is  little  credible  upwind  fallout  data  to 
serve  as  a  guide.  However,  a  correction  was  developed  that  is  reasonably 
consistent  with  data  available  for  shots  Koon,  Zuni  and  Bravo,  and  that 
provides  a  continuous  transition  to  the  low  yield  correction  at  W  =  10  KT. 

As  with  the  farfield  correction,  an  exponential  interpolation  in  log(ll) 
-  log(X)  space  is  used  to  compute  the  upwind-corrected  activity  deposited, 

D ( X)u ,  which  is  given  by 


D(X)u  =  D (X)exp 


bX 


r 

L 


1 .O-ex 


X<0 


where  D(X)  is  calculated  by  eqs.  (12a)  or  (12b),  and 


(IE 


b  =  0.0176 
c  =  570 

b  =  0 . 08045W~° ' 6G 
c  -  -8179.82  +  3800tnW 


W  <  10  KT 


W  >  10  KT 


(19 


^Though  one  would  guess  that  this  correction  should  be  a  function  of  wind 
speed,  there  is  not  enough  variation  of  v  among  the  cases  available  to 
allow  for  a  quantitative  evaluation  of  the  dependence. 


3.5  FALLOUT  TIME  OF  ARRIVAL 
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An  estimate  of  fallout  time  of  arrival  as  a  function  of  distance  from 
ground  zero  is  required  for  computation  of  maximum  effective  biological 
dose  and  turbulent  and  wind  shear  dispersion  of  the  nuclear  cloud  during 
atmospheric  transport.  By  time  of  arrival  we  mean  the  time  of  deposition 
of  the  first  fallout  at  distance  X  from  ground  zero  along  the  hotline. 

Time  of  arrival  ,  t  ,  is  estimated  by  means  of  the  following  simple 
model.  The  first  fallout  to  touch  ground  anywhere  does  so  at  onset  time 
t  (rq.  (1))  in  the  form  of  a  horizontally  distributed  parcel  centered  at 
X  -  vt().  We  take  the  radius  of  this  parcel  to  be  that  of  the  early  cloud, 
R^  (eq.  (2)).  Thus,  for  any  point  with  coordinate  X  <  vtQ  +  R.. ,  we  take 
t  =  tQ,  For  X  >  vtQ  +  R.j ,  we  take  ta  =  tQ  +  (X-vtQ-Ri )/ v.  The  geometry 
is  shown  in  Figure  6a. 

Figure  6b  shows  that  the  plot  of  t  vs  X  consists  of  two  straight 
lines  that  intersect  at  point  (tQ,  vtQ  +  R.).  We  desire  a  smooth  transi¬ 
tion  between  the  two  curves  rather  than  the  discontinuous  transition  shown. 
This  is  achieved  by  replacing  the  straight  lines  with  a  hyperbola  that  is 
asymptotic  to  both  lines  and  has  its  center  at  the  intersection  point  of 
the  lines.  Thus  t  (s)  is  calculated  from  the  equation 

a  — 


v  >  0.01  m  s"1  (20) 

Figure  7  shows  plots  of  tg  vs  X  for  an  effective  fallout,  wind  speed  of 
10  m  s  1  at  several  different  yields.  WSEG-10  results  also  are  shown.  Since 
according  to  WSEG-10  the  minimum  t  is  30  minutes,  we  see  that  for  low  yield 
shots  WSEG-10  grossly  overestimates  t  ,  and  hence  correspondingly  under¬ 
estimates  maximum  effective  biological  dose  (eq.  (32)). 

3.6  HORIZONTAL  SPREAD  OF  THE  NUCLEAR  CLOUD 

In  this  section  we  consider  horizontal  spread  of  the  nuclear  cloud 
before  we  account  For  dispersing  effects  of  atmospheric  turbulence  and  wind 
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Figure  6.  Basis  of  the  time-of-arri val  calculation 


fallout  time  OF  ARRIVAL,  to  (SI 


10  10  10  10  10  10 
downwind  distance  from  GROUND  ZERO,  X  (m) 


Figure  7.  Fallout  time  of  arrival  vs.  distance  from  ground  zero  for 
several  yields  as  computed  by  DNAF-1  and  WSEG-10. 
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shear.  Vertical  cloud  structure  is  included  implicitly.  We  have  shown 
elsewhere1'’7  that  close-in  fallout  patterns  from  surface  bursts  of  yields 
less  than  roughly  50  KT  are  dominated  by  fallout  from  the  cloud  stem.  Thus, 
for  a  low  yield  shot  we  cannot  make  the  conventional  assumption  that  close- 
in  fallout  cones  from  the  cloud  cap,  nor  that  it  begins  its  atmospheric 
transport  with  a  horizontal  spread  derived  from  the  stabilized  cloud  radius. 
Accordingly,  in  the  following  analysis  we  differentiate  between  stem  and 
cap  fallout. 

Three  critical  times  are  involved  here: 
t  fallout  onset  time  (eq.  (1)) 
t  fallout  arrival  time  (eq.  (20)) 

cl 

tg  time  of  ground  impact  of  a  nominal  particle  which 
begins  its  descent  at  the  stabilized  cloud  base,  zR. 

Time  tR  (s)  is  approximated  by  use  of  a  simple  relation  between  settling 
speed  of  water  drops  and  altitude  ,  which  for  this  purpose  is  found  to  apply 
well  enough  to  fallout  particles.  The  settling  speed  of  a  particle  at 
altitude  z,  f (z) ,  is 

f(z)  =  f(o)e^  (21) 

where  r,  =  2.90  X  10_J  m_1  and  from  Table  1,  f(o)  =  1  .6538  (m  s  1)  for  our 
nominal  particle.  Thus, 

r°  ( 

<■6  =  -  I  w(i-.7(tf|o|) 

where  Zq  is  given  by  eq0  (3). 

Horizontal  dimension  of  the  cloud  prior  to  atmospheric  transport  is 

specified  in  terms  of  the  standard  deviation  of  its  spread,  oc  (in).  Define 

a  yield  dependent  parameter,  o.  as 

w 
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ml 

m 

9 

■:r* 


°w 

=  Ri 

;  W  < 

10  KT 

(23a) 

°w 

=  Ri  +  (2,5Ri 

-  R-jMlogi 

0W  -  1 

)/ 2 

=  Ri(l  +  31 og 

i0W)/4;  10 

<  W  < 

1000  KT 

(23b) 

0 

w 

=  2.5  R. 

;  w  > 

1000 

KT, 

(23c ) 

where  R.  is 

given  by  eq. 

(2).  Then 

upwind 

and  in  the  region  of  ground 

zero 

we  have 


°c  -  aw 


;  X  <  vt  . 


For  fallout  from  the  stem  we  have 


=  “»  +  (t|^)(21  -  0  ;  X  >  vto  and  ‘a  '  l! 


and  for  fallout  from  the  cap  we  have 


°c  =  V2’  ta  -  tB 


(24) 


(25) 


(26) 


where  R  is  given  by  eq.  (5). 

Equations  (24)  and  (23a)  express  the  fact  that  onset  of  fallout  from  low 
yield  shots  is  early  enough  that  the  upwind  and  ground  zero  area  fallout  has 
essentially  the  spread  of  the  late  fireball.  Equations  (24)  and  (23c) 
account  for  the  fact  that  the  debris  from  high  yield  shots  is  carried  aloft 
rapidly,  which  causes  the  earliest  fallout  to  traverse  a  substantial  vertical 
path  and  thus  experience  substantial  horizontal  dispersion.  Equation  (23b) 
is  simply  a  linear  interpolation  in  log10(W)  between  eqs.  (23a)  and  (23c). 

Equation  (26)  sets  the  standard  deviation  of  horizontal  spread  of 
Fallout  in  the  stabilized  cloud  cap  at  one  half  of  the  stabilized  cloud 
radius,  and  eq.  (25)  provides  for  stem  fallout  via  a  linear  interpolation 
in  altitude  (in  terms  of  arrival  time)  between  the  base  and  top  of  the  stem. 
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3.7  TURBULENT  DISPERSION  OF  FALLOUT 


During  transport  from  its  initial  location  in  the  stabilized  cloud  to 
the  ground,  fallout  is  acted  upon  by  the  ambient  atmospheric  turbulence  such 
as  to  produce  additional  dispersion.  To  calculate  this  effect,  we  use  the 
scale  dependent  equations  of  Walton  which  require  specification  of  turbu¬ 
lence  level  in  terms  of  a  quantity  called  turbulent  energy  density  dissipation 
rate,  e.  Of  course,  e  will  depend  on  local  conditions  in  the  atmosphere, 

1  4 

but  Wilkins  has  found  that  c  can  be  approximated,  wit!)  surprisingly  con¬ 
sistent  accuracy,  by  a  simple  reciprocal  function  of  altitude.  Thus,  the 
variance  of  the  horizontal  spread  of  fallout  at  ground  level,  c/'  (m?),  not 
including  crosswind  dispersion  owing  to  wind  shear,  is  given  by 


2 

G 


1/3 .  \  3 

<f:>  "a  j 


109  m2 


(27a) 


•jlf 

o2  =  10C(3a2/3  +  2  <,->i/3ta  -  2000  )  ;  a?  >  1 0 V  (27b) 

where  t.  is  given  by  eq.  (20)  and  o,  is  calculated  as  described  in  the  pre- 

d  C 

ceding  section. 

Wilkins1  relation  for  e  is  v.  ^  0.03/z,  and  we  have  taken  for  our  average 
value,  <r.>  =  0„03/Zg.  Using  a  power  function  in  W  relation  for  which  is 
approximately  valid  over  the  entire  yield  range15,  we  obtain 

=  O.Ol6522W"0-10^  W  in  KT. 

?  9  2 

The  value  of  t  at  which  o  =  10  m  is  given  by 


*See  ref.  1,  sec.  3.3  for  a  more  complete  presentation  of  these  equations. 
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Thus  if  t  <  (t.  ).,  use  eq.  (27a);  otherwise  use  eq.  (27b). 

The  value  of  o  calculated  by  one  of  eqs.  (27)  is  used  in  one  of  eqs. 

(12)  to  calculate  crosswind  integrated  fraction  of  activity  deposited  at 
hotline  distance  X  from  ground  zero. 

At  this  point  in  the  presentation  we  have  discussed  how  to  determine 
all  of  the  quantities  needed  to  calculate  D(X)  via  one  of  eqs.  (12a)  or 
(12b).  Figure  8  shows  plots  of  D(X),  including  upwind  and  farfield  correc¬ 
tions,  at  every  other  decade  of  yield  for  an  effective  fallout  wind  of 
10  m  s“  .  WSEG-10  results  also  are  shown  for  comparison. 

The  most  obvious  difference  between  the  DNAF-1  and  WSEG-10  predictions 
is  the  much  sharper  peak  downwind  of  ground  zero  (GZ)  predicted  by  DNAF-1. 

For  the  higher  yields  this  peak  falls  off  more  rapidly  toward  G7,  according 
to  DNAF-1,  such  that  the  DNAF-1  GZ  activity  is  substantially  less  than 
predicted  by  WSEG-10,,  The  farfield  activity  curves  have  nearly  the  same 
shape,  as  expected,  though  they  are  significantly  displaced,  except  for 

4 

10  KT  for  which  case  they  are  essentially  coincident.  Upwind,  the  curves 
have  similar  shapes,  though  again  the  displacements  are  significant.  Shapes 
and  displacements  at  near  and  intermediate  downwind  distances  are  signifi¬ 
cantly  different. 

3.8  WIND  SHEAR  DISPERSION  AND  CROSSWIND  SPREAD  OF  THE  FALLOUT  PATTERN 
Following  in  principle,  but  not  in  detail,  the  procedure  of  Pugh  and 

If 

Galiano  ,  we  account  for  the  effect  of  vertical  wind  shear  on  crosswind 

2  2 

dispersion  variance  by  an  added  variance  increment,  os  (in'),  given  by 

as  =  [MzT  "  W10  ?  (28) 

Here  Sy  is  an  approximation  to  the  crosswind  component  of  vertical  wind 
shear,  determined  as  described  in  section  4.3,  and  the  other  quantities 
are  as  defined  by  eqs.  (3),  (4)  and  (20).  In  addition  to  being  a  very 
rough  approximation,  this  equation  is  somewhat  arbitrary  in  that  some  height 
difference  other  than  Zj  -  z^  could  have  been  used.  The  divisor  10  was 


T“n'rn  i  p 


"TITITT]" 


YIELD  =  100  KT 
WIND  =  10  m  s1 

DNAF-1  - 

WSEG-10 - 


DOWNWIND! 
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chosen  by  numerical  experimentation  to  give  good  comparisons  between  observed 
and  calculated  test  shot  fallout  patterns. 

Crosswind  (i.e. ,  Y  axis)  dispersion  of  the  fallout  pattern  is  provided 
by  multiplication  of  eqs.  (15)  and  (18)  by  a  Gaussian  function 


(29) 


where 


Oy?  -  O  +  Os  t  (30) 

2  2  2 

and  o  and  us  are  given  by  eqs.  (27)  and  (28)  (in  ). 

3.9  GAMMA  RAY  EXPOSURE  RATE  AND  MAXIMUM  EFFECTIVE  BIOLOGICAL  DOSE 

If  we  define,  as  usual,  X  to  be  the  distance  from  ground  zero  along 
the  direction  of  the  effective  fallout  wind  vector,  positive  in  the  down¬ 
wind  direction,  and  Y  to  be  perpendicular  distance  to  the  X  axis  in  the 
ground  plane,  then  the  11  f  1  hour  normalized  gamma  ray  exposure  rate 
(Roentgens  per  hour)  at  a  height  of  one  meter  above  a  point  X,  Y  on  the 
ground  is 


A(X,Y)  =  CKWfG(Y)D(X)  (31) 

where  C  is  a  scale  factor  (for  example,  see  Appendix  B),  Wp  is  fission 
yield  (KT),  K  =  6.9733  X  lo'1  (Roentgens  -  m?)/(hr  -  KT)**,  and  G(Y)  is 
given  by  eq.  (29).  D(X>0)^_p  is  given  by  eq.  (15),  while  D(X<0)  _u  is 
given  by  eq.  (18).  Fallout  maps  are  symmetrical  across  the  X  axis 
(i.e.,  A(X,Y)  =  A( X ,- Y) ) . 


*The  "normalized"  II  +  1  hour  exposure  rate  assumes  that  all  fallout  is  de¬ 
posited  at  H  +  1  hour,  regardless  of  whether  this  is  actually  the  case  or  not. 

**In  the  older,  more  familiar  units,  K  =  2692.4  (Roentgens  -  mi7)/ (hr  -  KT ) 


To  estimate  radiation  damage  to  people,  in  terms  of  short-term 

survi  vabi  1  i  ty  for  damage  assessment/ vulnerabi'I  ity  analysis  studies, 

a  quantity  here  called  "maximum  effective  biological  dose"  is  con- 
s 

ventionally  used  ,  This  quantity  is  designed  to  allow  for  effects  of 
a  continuing  exposure  of  ever  decreasing  intensity,  and  to  account  for 
some  coincident  repair  of  radiation  damage  by  the  human  body.  Follow¬ 
ing  a  theory  postulated  by  Blair16,  Davidson17  assumes  that  90  percent 
of  total  radiation  injury  is  reparable,  while  the  remaining  10  percent 
is  irreparable.  Further,  he  estimaies  that  for  humans  the  repair  rate 
is  about  0.1  percent  per  hour  of  the  residual  reparable  injury.  Taking 
fallout  gamma  radiation  exposure  rate  to  vary  with  time  according  to 
the  usual  t  1  “  approximation18,  Davidson  derives  an  equation  for  the 
ratio  of  biological  effective  dose  to  II  +  1  hour  exposure  rate  that  is 
a  function  of  two  variables:  time  of  arrival  of  fallout  (or  time  of  entry 
into  the  fallout  field),  and  time  of  exit  from  the  fallout  field.  This 
equation  has  been  evaluated  numerically,  and  when  plotted  against  exit 
time  for  specified  t, ,  the  curve  is  found  to  have  a  maximum:  the  late- 

a 

time  falloff  in  effective  biological  dose  being  caused  by  combined 
effects  of  damage  repair  and  deray  of  exposure  rate  intensity.  The 
maximum  in  this  curve  gives  the  quantity  called  maximum  effective 
biological  dose,  M(X,Y),  and  if  we  assume  that  residence  in  che  fallout 
field  is  from  tQ  to  at  least  the  time  of  the  maximum,  it  is  a  function  only 
of  A(X,Y)  and  tfl ,  where  A ( X ,  Y )  is  a  simple  multiplier. 

The  numerical  calculations  necessary  to  define  the  ratio  M ( X  ,  Y ) / 

A(X,Y)  as  a  function  of  tg  have  been  done  by  the  DoD  Command  and  Control 
Technical  Center  and  simple  functions  have  been  fitted  to  the  results  to 
give  the  following  "quick  approximation"  equations: :> 

M(X,Y)  =  A(X,Y)(a0  +  a^  +  a^Z) 


4b 


(32) 


where 


B  =  0.8685833  ln(tj 

a 

a  -  1  5.2891 
o 

a-j  =  -2. 903225 
a2  =  0.1662315 


t  <  1157.9  s 

d 


B  =  2  1 n(t  ) 

a 

a0  =  4.6182 
a,  =  -0.53587 
a^  =  0.016923 


t  >  1157.9  s. 

a  — 


m : 


■  ,i 

■  1 


As  discussed  in  sections  3.1  arid  3.2,  a  Gaussian  dispersion  function 
for  deposited  fallout  is  preferred  for  both  the  alongwind  and  crosswind  direc¬ 
tions.  In  this  model  a  Gaussian  crosswind  function,  G(Y),  is  used,  but  the 
alongwind  function,  F(X,t)  (eq.  (10)),  is  non-Gaussian .  These  functions  are 
compared  in  Figure  4.  A  consequence  of  this  inconsistency  is  that  the  fall¬ 
out  pattern  is  always  asymmetric,  even  for  zero  wind,  in  which  case  all 
activity  contours  should  be  circles  centered  at  ground  zero.  Specifically 
for  v  =  0,  eq.  (31)  becomes,  if  we  omit  the  upwind  and  farfield  corrections, 


A(X,Y);.,0 


KWr 


exp 


Ml 


o?y[2?[]  *  (£)2] 


,3 


which  obviously  cannot  give  circular  contours  for  A  =  constant. 

In  practical  terms  this  defect  in  the  model  is  of  little  consequence. 
This  is  because  a  zero  effective  fallout  wind  is  physically  unacceptable. 
Indeed,  it  has  been  iround  that  for  other  reasons  (see  sec.  4.2),  the  mini¬ 
mum  acceptable  value  of  v  is  about  0.5  in  s"  . 


lA 
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4.  USE  OF  WIND  DATA 


4.1  GENERAL  CONSIDERATIONS 

Wind  data  are  used  to  determine  two  essential  model  parameters:  the 
effective  fallout  wind  vector,  v,  and  the  crosswind  shear  parameter,  Sy. 
Use  of  the  magnitude  of  the  effective  fallout  wind  vector,  v,  is  described 
throughout  section  3,  and  use  of  Sy  is  explained  in  section  3.8. 

The  code  accepts  wind  data  in  two  forms:  either  the  user  can  specify 
v  and  Sy  directly  or  he  can  supply  a  single  vertical  profile  of  wind  vec¬ 
tor  data,  in  which  case,  the  code  computes  v  and  Sy  from  these  data.  In 
this  chapter  we  describe  these  computations. 


4.2  EFFECTIVE  FALLOUT  WIND 

The  code  accepts  a  single  vertical  profile  of  wind  vectors,  each  vector 
representing  the  wind  speed  and  direction  at  a  specified  altitude.  As  is 
described  in  detail  in  section  6.3,  considerable  flexibility  is  allowed  in 
terms  of  form  and  format  of  the  input  data. 

After  some  preprocessing  (subroutine  INWIND),  the  data  are  stored  in 
tabular  form.  There  are  four  tables  which  contain  the  following  data: 
zi’  ^E.i’  UN  i  and  zb  i*  Here  zi  1S  the  altitude  (in  above  ground)  at  which 
wind  vector  components  IJ^  ^ ..  arc  defined*,  i  is  the  table  entry  (i.e., 
wind  stratum)  index  (i  =  1,  2  -  -  -  I),  and  z^  -j  is  the  base  altitude  (in  above 
ground)  of  the  i^h  wind  stratum  defined  as 


2  ^zi-l 


*Note  that  it  is  standard  practice  to  measure  surface  wind  at  an  elevation 
of  1 0  meters . 


47 


with  z,  ,=0.  Up  .  is  the  wind  component  along  the  west-east  axis,  positive 
toward  the  east,  and  .  is  the  wind  component  along  the  south-north  axis, 
positive  toward  the  north  (m  s  1). 

Strictly  speaking,  altitude  should  be  relative  to  mean  sea  level  (MSI). 
However,  in  most  cases  MSL  can  be  replaced  by  ground  level  (GL)  without  sub¬ 
stantial  error,  and  in  practice  this  substitution  will  be  implicit,  in  most 
land  surface  burst  predictions,  as  it  is  in  the  cases  of  the  predictions  of 
the  Nevada  Test  Site  shots  discussed  in  section  5.  The  code  provides  for 
adjustment  of  altitudes  to  be  relative  to  GL  even  though  they  may  be  input 
relative  to  some  other  origin. 

Effective  fallout  wind  is  a  weighted-average  wind,  the  average  being 
taken  between  the  stabilized  cloud  cap  center  height  and  the  surface,  and 
the  weighting  being  taken  according  to  settling  time  of  the  nominal  par¬ 
ticle  (sec.  2.4)  through  each  wind  stratum.  The  calculations  are  done  in 
subroutine  EFWIND. 

Define  U.  to  be  the  wind  vector  in  the  itl1  stratum.  Then  the  effective 
fallout  wind  is 


J-l 

£{“<(  *!>,«♦>  -zb,l) 

i/'Up  ♦  ii2'i 
J 

(zc  '  zb,J 

)/f(zp 

J-l 

^  (  zb, i+1  zb 

,()/f(zp  ♦  ( 

zc  ‘  zb,j) 

/f(zj) 

where  the  summation  begins  at  the  ground,  z-  and  z^  •  are  as  defined  above, 
but 


zc  “  (ZT  +  zl))/2 

zj  ■  (zc  '  zb,j)/2  ;  ZM  -  Zc 
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and  U  1  is  the  wind  vector  at  altitude  z '  as  determined  by  linear  interpolation. 
z  j  J 

f(z)  is  the  settling  speed  of  the  nominal  particle  at  altitude  z  determined 
by  linear  interpolation  in  Table  1. 

Actually,  the  code  uses  the  magnitude  of  the  effective  fallout  wind, 


(35) 


and  the  sine  and  cosine  of  its  direction  angle  $  defined  as 


sin  ij>  =  VjVv 
cos  <J>  =  vN/v 


(36) 


where  v^  and  v^  are  the  easterly  and  northerly  directed  components  of  v. 

Theoretical  and  practical  considerations  impose  a  lower  limit  on  the 
acceptable  value  of  v„  While  occasionally  a  calm  condition  may  bo  observed 
at  the  surface,  this  is  never  the  case  throughout  the  transport  air  space, 
and  therefore  a  zero  value  for  v  is  never  acceptable.  Very  low  values  of 
v  may  cause  certain  unrealistic  results  to  appear:  for  example,  the  upwind 
hotline  activity  may  fall  off  less  rapidly  than  the  downwind  activity.* 

Accordingly,  the  code  will  not  accept  a  value  of  v  less  than  0.5  m  s 
An  input  value  of  v  =  0  is  used  as  a  flag  to  signal  input  of  a  vertical 
profile  of  wind  data.  When  the  code  encounters  a  value  v  loss  than  0.5  in  s 
(which  is  not  interpreted  as  the  value  0.0  used  to  signal  input,  of  the  ver¬ 
tical  profile)  this  value  is  printed  along  with  a  comment,  and  v  is  reset 
to  0.5  m  s' 1 . 


vf[rfs  anemia  lous  behavior  is  caused  by  interaction  of  several  features  of 
the  code,  lirsl,  the  horizontal  variance  of  deposited  fallout,  a  (secs. 
3.6  and  3.7)  is  held  constant  upwind  of  ground  zero,  whereas  it  increases 
downwind.  Second,  the  upwind  correction  (sec.  3.4)  is,  unfortunately,  not 
a  function  of  v,  but  was  determined  from  test  shot  results  for  which  v  is 
always  substantially  greater  than  zero.  Consequently,  both  of  these  fea¬ 
tures  depend  on  use  of  realistically  large  values  of  v  to  give  realistic 
results. 
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4.3  SHEAR  PARAMETER 


Vortical  wind  shear  is  defined  .is 


'  dU 
b  "  d  t 


where  II  is  wind  and  z  is  the  vertical  coordinate.  We  make  the  customary 
assumption  that  advective  transport,  will  overwhelm  effects  of  shear  dis¬ 
persion  in  the  alongwind  direction,  and,  iheref'nre.  we  are  interested 
only  in  the  crosswind  component  of  S,  Sy . 

In  this  model  ,  Sy  is  taken  to  be  the  root-mean- square  value  of  the 
crosswind  components  of  aLI/az  computed  at  intervals  of  Az  =  (zy  -  z.?)  from 
the  cloud  ton  to  the  ground.  The  final  Az  value  is  adjusted  as  required 
to  avoid  reaching  below  the  ground. 

In  terms  of  the  variables  defined  in  the  preceding  section. 


irE  H(V, 


c°s*  (uEj^  -  II, 


z  .  -  z  • 


1/p 

i  )  (38) 


Here  the  summation  begins  at  the  cloud  top  such  I  ha  I  we  have 


Zj  '  ZT 


-  (j  -  1  )(zT  -  zH  ]  ;  j  1,2, - ,  K-l 


zK  -  0. 


5.  VALIDATION 


5.1  DISCUSSION  OF  RESULTS 

Predictions  are  compared  with  observed  H  +  1  hour  normalized*  expo¬ 
sure  rate  maps  for  the  first  five  test  shots  described  in  Table  2.  For 
the  sixth  shot,  Bravo,  there  are  not  enough  observed  data  to  construct  a 
complete  fallout  map.  Thus,  for  this  case  we  compare  our  prediction  against 
a  special  "reconstruction"  calculation  made  by  the  Naval  Radiological  Defense 
Laboratory  shortly  after  the  event19. 

Three  methods  of  comparison  of  fallout  patterns  are  used: 

1.  Visual  comparison  of  contour  maps. 

2.  Comparison  of  contour  areas,  and  hotline  lengths  and  azimuths.** 

TABLE  2 

TEST  SHOT  DATA 

Total  Fission  Altitude 


Yield 

Yield 

HOB 

of  GZ 

Shot 

(KT) 

(KT) 

("0 

(m) 

Site 

Johnie  Boy 

0.5 

0.5 

-0.584 

1570.6 

NTS+ 

Jangle-S 

1.2 

1.2 

1.067 

1284.7 

NTS 

Small  Boy 

1  ow 

- 

3.048 

938.2 

NTS 

Koon 

150. 

- 

4.145 

0.0 

Bi kini 

Zuni 

3380. 

- 

2.743 

0.0 

Bi kini 

Bravo 

15000. 

- 

2.134 

0.0 

Bi  ki ni 

•f 

Nevada  Test  Site 


A  "normalized"  exposure  rate  map  is  constructed  on  the  assumption  that  all  local 
fallout  is  down  at  the  specified  time,  regardless  of  its  actual  deposition  time. 

Hotline  length  is  defined  as  the  furthest  distance  from  ground  zero  on  a  contour, 
and  hotline  azimuth  is  the  angle,  measured  clockwise  from  north,  to  the  point  of 
furthest  distance  from  ground  zero  on  a  contour. 
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2  0 

3.  The  Row! and-Thompson  Figure-of-Merit  (FM)  which  is  a  measure  of 
contour  overlap.  (See  Appendix  C.) 

These  are  roughly  in  order  of  importance. 

Statistical  data  are  in  Table  3  and  the  contour  plots  are  on  pp.  58 
through  80.  Contours  were  drawn  by  a  30-inch  Calcomp  plotter,  and  each 
observed-predicted  pair  are  to  the  same  scale.  Contour  maps  and  statisti¬ 
cal  data  are  included  for  predictions  by  DFLFIC  and  WSEG-10  as  well  as  by 
DNAF-1. 

Prediction  accuracy  is  seen  to  be  good.  Perhaps  the  best  quantita¬ 
tive  measure  of  accuracy  is  provided  by  the  mean  absolute  percent  error, 

E,  which  for  n  observed-predicted  data  pairs  is 


E  *  ™  Ehbs.i  -  * 


/  X 

pred,i  obs,i 


Values  of  E  for  each  prediction  (excluding  Bravo)  by  each  of  the  three 
models  are  given  under  the  solid  lines  in  Table  3.  The  values  in  paren¬ 
theses  are  computed  with  the  data  for  the  highest  level  contours  excluded. 
The  highest  level  contours  are  particularly  difficult  to  predict  since 
usually  they  are  dominated  by  the  region  most  affected  by  induced  activity 
in  the  ground  and  throwout  from  the  crater,  neither  of  which  are  addressed 
by  the  fallout  models.  Overall  mean  absolute  percent  errors  are  given  in 
Table  4.  (Bravo  prediction  data  are  excluded.)  As  one  might  expect, 

DNAF-1  errors  are  'intermediate  between  those  of  DFLFIC,  which  are  best, 
and  WSEG-10,  which  are  worst,  though  the  differences  between  the  DNAF-1  and 
DFLFIC  errors  are  less  than  between  DNAF-1  and  WSEG-10.  Note  that  the  most 
obvious  problem  with  the  WSEG-10  predictions  is  a  tendency  to  overpredict 
the  low  level  contours  at  the  expense  of  the  higher  levels,  to  the  extent 
that  frequently  the  higher  level  contours  are  completely  absent. 
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OVERALL  ML' AN  ABSOLUTL  PERCENT  ERRORS* 


Contour  Area 

DNAF-1  66(08) 

DELEIC  62(4?) 

WSEfi- 10  117(90) 

Thu  ricjuro-of -Merit  (EM)  results 
capabilities  for  the  models.  This  is 
and  we  have  concluded  that  in  its  pres 
useful  measure  of  prediction  capabilit 


Hot! ine  Length 

43(35) 

3? (?6) 

51(45) 

do  not  show  a  consistent  order  of 
typical  of  past  experience  as  well, 
out  form,  EM  does  not  provide  a  very 
y.  Details  are  given  in  Appendix  C. 


5.2  DISCUSSION  01  THE  1  ESI  SHOT  DATA  AND  PREDICTIONS 

The  three  low  yield  shots  were  executed  at  the  Nevada  Test  Site,  and 
their  fallout  patterns  were  measured  over  land,  for  this  reason,  observi  d 
patterns  for  these  shots,  though  not  highly  accurate,  may  be  considered  to 
be  superior  to  the  patterns  of  the  high  yield  shots  which  were  executed  on 
bikini  Atoll  in  the  South  Pacific.  Not  only  are  the  fallout  fields  of  the 
high  yield  shots  very  large,  which  adds  to  measurement  problems,  but  most 
of  the  fallout  from  these  shots  fell  into  water.  Even  so,  most  of  the  Koon 
pattern  area  was  covered  by  an  array  of  fallout  collection  stations,  so  this 
pattern  is  probably  reasonably  accurate,  Zuni,  on  the  other  hand,  is  a 
special  case.  The  fallout  pattern  list'd  here  is  exclusively  downwind  of  the 
atoll  and  was  determined  by  an  oceanographic  survey  iih  i  hod  that  was  known 
to  be  inaccurate.  The  close-in  pattern  in  the  region  of  the  atoll  is 
available,  but  contains  no  closed  contours  so  '  .  could  not  be  used  hero;  thus 
the  it  igh  acti v  i  ty  portion  of  the  observed  pattern  for  this  shot  is  ignored 


*Valuer>  in  parentheses  are  calculated  with  da  La  for  the  highest  level 
contours  excluded. 
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and  this  alone  must  account  for  a  substantial  portion  of  the  disagreement 
between  observation  and  prediction  for  this  shot,  particularly  with  regard 
to  contour  areas  and  contour  overlap  (Table  3).  As  already  mentioned,  we 
have  no  observed  pattern  for  the  Bravo  shot.  In  addition,  we  have  the 
following  problem. 

DNAF- 1  and  DELFIC  predictions  for  the  high  yield  shots  are  expected 
to  be  inferior  to  those  for  the  low  yield  shots.  This  is  because  the  high 
yield  shots  were  detonated  over  coral  soil,  and  in  the  cases  of  Zuni  and 
Bravo,  large  but  uncertain  amounts  of  sea  water  were  lifted  by  the  clouds. 

The  particle  size  distribution  used  for  these  predictions  is  typical  of 
fallout  produced  from  the  siliceous  soil  found  at  the  Nevada  Test  Site, 

We  have  not  succeeded  in  developing  a  distribution  appropriate  for  coral 
and  coral-sea  water  mixtures. 

DNAF-1  predictions  were  made  using  the  H  hour  winds  tabulated  in  refer¬ 
ence  7.  DELFIC  predictions  were  made  using  all  of  the  reference  7  wind 
profiles,  from  H  hour  onward  in  time.  WSEG-10  calculations  were  done  using 
v  and  Sy  values  supplied  by  the  DoD  Command  and  Control  Technical  Center  as 
determined  by  them  from  the  H  hour  wind  profiles;  these  data  also  are  tabula¬ 
ted  in  referen  7  (Appendix  A.3).  For  shots  Small  Boy  and  Bravo,  the 
published  wind  data  have  been  found  to  be  not  pertinent  to  transport  of  the 
nuclear  clouds.  For  both  of  these  cases,  we  have  used  reconstructed  wind 
data:  for  Small  Boy  the  reconstruction  is  described  in  Appendix  B  of 
reference  7,  ..nd  for  Bravo  we  have  used  the  winds  developed  by  Dean  and 
Olmstead.  Values  of  v,  f  and  Sy  computed  for  the  DNAF-1  predictions  are 
given  in  Table  5. 


TABLE  5 
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EFFECTIVE  FALLOUT  WINDS  AND  SHEAR  PARAMETERS 
COMPUTED  FROM  H  HOUR  WIND  PROFILES  FOR  USE  BY  DNAF-1 


Test 

,  V-2, 

(ill  S  ) 

<i> 

i£ll 

Shot 

(deg. ) 

Johnie  Boy 

6.0 

-  8.6 

0.00323 

dangl e-S 

13.1 

14.6 

0.00311 

Sma 1 1  Boy 

3.8 

64.0 

0.00066 

Koon 

6.2 

11.3 

0.00133 

Zuni 

4.9 

-20.0 

0.00225 

Bravo 

5.8 

93.6 

0.00044 

5„  3  OBSERVED  AND  PREDICTED  FALLOUT  PATTERNS 

Contours  are  in  units  of  Roentgens  per  hour  for  gamma  radiation  at  a 
height  of  one  meter  above  ground  at  H  +  1  hour.  All  activity  is  assumed 
to  be  deposited  at  H  +  1  hour.  For  all  but  the  Zuni  shot,  for  which  fall¬ 
out  activity  was  measured  by  an  oceanographic  method,  predicted  activities 
are  multiplied  by  a  combined  ground  roughness-instrument  response  correction 
factor  of  0.5.  (See  Appendix  B.) 

Observed  and  predicted  patterns  for  each  case  are  plotted  to  the  same 
scale.  North  is  up  the  pages  and  east  is  across  the  pages  toward  the  right. 
Visual  comparisons  are  best  made  by  superimposing  electrostatic  copies  of 
the  plots. 
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6.  COMPUTER  CODE 


■  f 


6.1  GENERAL  DISCUSSION 

Coding  is  entirely  in  ANSI  FORTRAN  and  complete  listings  are  presented 
in  Appendix  D.  Glossaries  of  mnemonics  are  in  subroutines  DNAF  and  INWIND, 
and  the  Appendix  A  glossary  relates  symbols  used  in  the  report  text  to  the 
FORTRAN  mnemonics.  The  code  was  originally  developed  and  exercised  on  the 
CDC  6600  computer,  but  also  has  been  extensively  tested  arid  exercised  on 
the  Honeywell  Information  Systems  6080  computer  used  by  the  DoD  Command 
and  Control  Technical  Center.* 

Data  is  input  via  system  unit  5  and  output  via  unit  6.  These  unit 
specifications  can  be  changed  by  changing  the  values  of  parameters  IN  and 
10  (lines  69  and  70,  subroutine  DNAF).  No  other  peripheral  storage  units 
are  addressed. 

All  subroutine  and  function  names  are  listed,  with  a  brief  description, 
in  Table  6. 

There  are  two  optional  modes  of  code  operation: 

Mode  1.  H  +  1  hour  exposure  rate  map  preparation 

Mode  2.  H  +  1  hour  exposure  rate,  and  maximum  effective 
biological  dose  at  user  specified  points. 

Map  calculation  is  discussed  in  detail  in  section  6.2.  For  the 
Mode  2  calculation,  the  user  simply  specifies  ground  points  in  terms 
of  X,Y  coordinates  relative  to  ground  zero,  and  the  H  +  1  hour  normalized 
exposure  rates  (Roentgens  nr  )  and  maximum  effective  biological  doses 
(Roentgens),  at  a  height  of  one  meter  above  the  point  -,e  printed. 

Wind  data  may  be  specified  in  cither  of  two  ways: 

1.  Fffect;,e  fallout  wind  sneed,  v  (m  s  *),  wind  direction, 

rf>  (donees  clockwise  from  north  in  the  direction  of  the  wind 
vector),  and  ci osswind  shear  parameter,  Sy  (s-1),  are  input 
directly. 


*In  acuition  to  the  FORTRAN  code,  the  Texas  Instruments  TI  59,  a  sophisti¬ 
cated  pocket  calculator,  has  been  programmed  to  compute  A(X,0)  (see  eq.  31)). 
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TAHI  F  6 


i>;  ;■ ciui’TiON  or  qnaf-i  com  suiimmi  ini  s  and  functions 


Subrout i no 
or  Function 

!)NAF 

SCI.  00 

INWINll 

IT  WIND 
nVWNI) 


ONSLT 

DOST 

I1NAT  I 

SI  OTA 


Descri ption 

Executive.  Controls  calculation  flow  and  calls  most  of  the  followin']  urn!]rains. 

Computes  stabilized  cloud  heights  and  radius  ulus  fireball  radius: 

z(r  zT,  NCj.  Rj. 

Heads  in  and  processes  a  vertical  profile  of  wind  data  if  the  effective 
fallout  wind  speed,  v,  is  not  specified  by  input. 

Computes  effective  fallout  wind,  v.  from  data  input  via  subroutine  INWINIl. 

Computes  the  wind  shear  parameter ,  Sy,  from  the  wind  data  input  via  Suhrout  ini' 
INWIND. 

Computes  fallout  onset  time,  t  . 

Computes  the  factor  which  when  i-ultiplied  by  |i  i  1  hour  exposure  rate  <]  i  ves 
the  maximum  effective  biological  dose,  M(X,Y). 

Computes  erosswind-integrated  H  *■  I  hour  activity  fraction  deposited  at  distance 
X  aioncj  the  hotline  from  ground  zero. 

Compute',  deposition  variance,  ■  .  crosswind  pattern  standard  deviation,  ,iy, 
and  fallout  arrival  time,  t  .  all  a  i.  a  Innnw  ind  distance  X  from  ground  zero. 


MAIN1 

SITMT 


T  HI’I 


Prints  a  two-dimensional  map  of  h  '  1  hour  exposure  rate. 

Sets  iallout.  map  boundaries  and  grid  intervals  for  the  user.  This  subroutine 
is  talk'd  only  if  parameter  DSI  TMP  is  false,  and  tile  resulting  map  is  intended 
tu  provide  the  user  with  a  preliminary  look  at  the  fallout  pattern. 

Linear  interpolation  between  table  entries. 


EkllOlf 


Error  condition  printout. 


nx 


2.  A  vertical  profile  of  winds  are  input  from  which  the  code  com¬ 
putes  v  and  Sy  as  described  in  sections  4.2  and  4.3. 

The  effective  fallout  wind  speed  is  limited  to  values  greater  than  or 
equal  to  0.5  m  s-1.  A  zero  input  value  is  used  to  signal  input  of  a  verti¬ 
cal  profile  of  wind  data.  If  a  value  of  v  <  0.5  m  s’1,  other  than  the  signal 
value  of  0.0,  is  encountered,  a  comment  is  printed  and  v  is  reset  to  0.5 
in  s”1.  (See  the  end  of  section  4.2.) 

6.2  FALLOUT  MAPS 

Selection  of  the  mode  1  calculation  option  (logical  parameter  IFMAP, 
card  3,  is  true)  causes  calculation  of  II  +  1  hour  normalized  exposure  rates 
(Roentgens  hr  1 )  at  a  height  of  one  meter  above  a  two-dimensional  array  of 
points  on  the  ground.  The  points  are  spaced  at  grid  intervals  DGX  and  DGY 
(meters)  in  the  x  and  y  coordinate  directions.  Here  the  x  axis  is  in  the 
west- to-east  direction  (positive  east,  of  ground  zero),  and  the  y  axis  is 
in  the  south-to-north  direction  (positive  north  of  ground  zero). 

Ordinate  values  are  output  in  ^octangular  arrays,  arid  it  is  assumed 
that  the  arrays  arc  printed  by  a  standard  line  printer.  The  x  axis  is 
across  the  printed  page,  west-to-east  from  left-to-right,  and  the  y  axis 
is  up  the  page,  south- to-north  from  bottom-to-top.  A  two-line,  power  of 
ten  format  is  used.  Thus,  the  activity  for  each  point  appears  as 

+NNNNN 

V.VVV 

which  is  interpreted  as 

iNNfoNN  i 

V.VVV  X  10-  (Roentgens  hr’  ). 

Printing  is  done  in  units  of  map  strips,  each  strip  consisting  of  a  suffi¬ 
cient  number  of  connected  printer  sheets  to  cover  the  entire  v  axis  range. 


Each  row  (across  the  page)  of  each  strip  contains  a  maximum  of  nineteen 
ordinate  points  in  the  x  direction.  Enough  strips  are  produced  to  cover  the 
entire  x  axis  range.  These  ..trips  can  be  attached  side-by-side  such  as  to 
construct  the  complete  map,  and  contours  may  then  be  drawn  by  hand,  x  and 
y  values  are  printed  at  regular  intervals  on  each  strip. 

If  logical  parameter  USETMP  is  true  (card  3),  the  user  must  specify 
boundaries  and  grid  increments  for  the  map.  Otherwise,  the  code  sets  these 
parameters  (via  subroutine  SETMP)  automatically,  using  both  yield  and  wind 
as  criteria.  The  result  is  a  small,  rather  poorly  resolved  map  which  may 
not  satisfy  the  particular  needs  of  the  user.  It  is  intended  to  provide 
the  user  with  a  preliminary  view  of  the  map.  From  the  information  gained 
from  this  quick  look,  the  user  may  devise  his  own  map  specifications  as 
described  next. 

To  specify  his  own  map,  the  user  must  supply  the  following  information: 

1.  Logical  parameter  USETMP  -  -TRUE*  (card  3,  sec.  6.5) 

2.  Coordinates  (x^,  Ym^n  O'”  XMIN,  YM IN)  of  the  southwest  corner 
of  the  map,  and  (x  .  ym, u  or  XMAX,  YMAX)  of  the  northeast 
corner  of  the  map.  (card  6) 

3.  Grid  increments  DGX,  or  DGX  and  DGY ,  in  the  x  and  y  axis  directions 
(card  6) 

If  only  DGX  is  specified,  the  code  computes  DGY  such  as  to  produce  a 
spatially  undistorted  map  on  a  standard  line  printer:  that  is,  one  with  10 
characters  per  inch  across  the  page  and  6  characters  per  inch  down  the  page. 
To  adjust  for  nonstandard  character  spacing,  parameters  IH  and  IV  (tines 
69  and  70  in  subroutine  DNAF  and  13  and  14  in  subroutine  SETMP)  must  be 
changed. 

The  values  of  XMIN  and  YMIN  specified  on  card  6  should  be  one  grid 
increment  less  than  the  values  actually  expected  on  the  printed  map. 


The  code  presented  in  Appendix  D  provides  for  a  maximum  of  5000  map 
points.  If  the  map  specified  by  card  6  input  requires  more  than  this 
number  of  points,  DGX  and  DGV  are  adjusted  such  that  no  more  than  5000 
map  points  are  required,  a  comment  to  the  effect  that  the  adjustment  has 
been  made  is  printed,  and  calculation  then  proceeds.  Map  point  ordinates 
are  stored  in  array  0MAP,  and  parameter  NMAP  is  used  as  a  variable  dimen¬ 
sion  for  OMAP.  To  change  the  maximum  number  of  points,  change  the  dimension 
of  the  OMAP  array  and  the  value  of  NMAP  as  desired  (lines  67,  69  and  70  of 
subroutine  DiNAF). 

6.'  INPUT  OF  WIND  PROFILE  DATA 

If  parameter  WIND  (which  is  the  effective  fallout  wind  speed,  v)  on 
card  2  is  zero,  the  code  calls  subroutine  INWIND  which  reads  in  a  vertical 
profile  of  wind  data  as  described  by  cards  4a-4n  (sec.  6.5).  This  input 
is  designed  for  maximum  versatility. 

Card  4a  is  used  to  specify  whether  the  wind  data  are  input  in  resolved 
form  (i.e.,  vector  components  in  the  easterly  and  northerly  directions), 
or  in  terms  of  speed  and  direction  angle. 

Card  4b  is  an  object-time  format  to  be  used  to  read  the  data. 

Card  4c  contains  scale  and  translation  factors,  and  card  4d  contains 
data  field  pointers.  Cards  4e  to  4n-l  each  contain  altitude  and  wind  vec¬ 
tor  data  for  a  wind  stratum,  and  the  last  card,  4n,  is  the  data  set  termina¬ 
tor. 

Use  of  the  scale-translation  data  and  the  field  pointers  deserves  some 
explanation. 

In  Cuinb  i na L  ion  with  trie  object-time  format,  the  field  pointers,  Nl, 

N2,  N3,  allow  any  arrangement  of  the  data  on  the  cards;  the  only  restriction 
being  that  all  of  the  data  cards  have  the  same  arrangement.  Each  card  con¬ 
tains  three  items  of  data; 
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1.  Altitude  of  the  wind  measurement 

2.  x  component  of  the  wind  vector 

3.  y  component  of  the  wind  vector 


FORM  =  RESO 


2.  Wind  direction  angle 

3.  Wind  speed 


FORM  =  METE 


Pointer  N1  is  associated  with  the  altitude,  N2  with  the  x  wind  component 
or  direction  angle,  and  N3  with  the  y  wind  component  or  wind  speed.  Collec¬ 
tively  Nl,  N2,  N3  consist  of  some  permutation  of  the  integers  1,  2,  3.  The 
object-time  format  specifies  three  data  fields,  and  Nl  specifies  which  of 
the  three  contains  1  he  altitude,  N2  specifies  which  contains  the  x  wind 
component  or  the  direction  angle,  etc.  For  example,  if  we  have  Nl  =  1, 

N2  =  3,  N3  =  2,  then  the  altitude  is  in  the  first  field  (from  the  left), 
the  x  wind  component  or  direction  angle  is  in  the  third  field,  and  the  y 
wind  component  or  wind  speed  is  in  the  second  field. 

The  scale  and  translation  data  input  via  card  4c  allows  input  of  data 
in  any  units,  and  allows  certain  common  data  translations  to  be  made.  Scale 
factors  are  in  fields  1-3  and  translations  are  in  fields  4-5  of  card  4c. 
After  application  of  these  translations  and  scale  factors,  altitude  must  be 
in  meters  with  origin  at  the  ground,  and  the  wind  must  be  expressed  in  terms 
of  components  (m  s-1)  in  the  x  (easterly)  and  y  (northerly)  directions.  If 
it  is  used,  the  wind  direction  angle  must  be,  after  scaling  and  translating, 
the  angle  of  the  wind  vector  measured  clockwise  from  north. 

Scale  and  translation  data  are  road  from  card  4c  into  array  SCALE(  ), 
The  three  data  items  on  each  wind  data  card  are  read  into  array  AP(  ).  For 
the  altitude,  SCALE(4)  is  a  translation,  which  is  applied  before  scaling,  to 
adjust  the  origin  to  be  at  ground  level,  and  SCALF(l)  is  a  scale  factoi 
used  to  adjust  tne  units  to  meters.  Thus  the  altitude,  ZCH,  is 


ZCH  =  (AP(N1)  +  SCALE (4)  )*SCALF.  ( 1 ) . 

If  FORM  =  RESO  (card  4a),  the  wind  vector  data  are  input  in  component, 
form.  In  this  case,  the  only  other  card  4c  datum  used  is  SCALE(2),  which 
is  a  scale  factor  applied  to  both  components  to  adjust  units  to  meters  per 
second.  Specifically,  we  have  for  the  x  and  y  wind  components,  WX  and  WY, 

WX  =  AP(N2)*SCALE(2) 

WY  =  AP(N3)*SCALL'(2). 

If  FORM  =  METE  (card  4a),  the  wind  vector  data  are  input  in  terms  of 
direction  angle  and  speed.  As  with  the  previous  case,  SCALE ( 2 )  is  used 
to  scale  the  wind  speed  to  units  of  meters  per  second.  SCALE (3)  is  a  scale 
factor  used  to  convert  the  direction  angle  to  degrees.  $CALE(5)  is  an  angle 
translation  (i.e.,  rotation)  which  is  input  in  the  same  units  as  is  the 
direction  angle.  The  angle  translation  calculation  is  set  up  to  convert 
the  angle  from  the  conventional  meteorological  specification  of  direction 
from  which  the  wind  is  blowing,  to  direction  toward  which  the  wind  is 
blowing;  in  other  words,  the  code  automatically  rotates  the  wind  angle 
through  180°  unless  this  is  circumvented  by  appropriate  specification  of 
SCALE(b).  Specifically  we  have 

WX  AP(N3)*SCALE(2)*SIN(-1-Jq(/\P(N2)*SCALE(3)  +  TRNS) ) 

WY  =  AP(N3)*SCALE(2)*C0S(ygg(AP(N2)*SCALi:(3)  +  TRNS) ) 

where 

TRNS  -  SCAL L ( b ) *5CAI. E£ ( 3 )  -  180. 


Default  values  for  the  scale  factors,  SCALE(l)  through  SCALE(3),  are 
unity.  That  is,  if  any  or  all  of  the  first  three  fields  on  card  4c  are  left 
blank,  the  code  sets  the  corresponding  scale  factors  to  unity. 

The  wind  data  should  extend  at  least  up  to  the  cloud  top  height,  Zy 
(eq.  (4)),  but  if  it  does  not,  the  code  simply  assumes  that  the  data  in  the 
highest  wind  stratum  is  constant  up  to  Zy. 

The  wind  data  are  sorted  by  the  code  and  arranged  in  sequence  of  in¬ 
creasing  altitude;  thus  it  is  not  necessary  that  they  be  input  in  any  par¬ 
ticular  order. 

An  example  of  a  common  situation  will  illustrate  how  to  set  up  the 
data  cards.  Suppose  the  wind  data  are  available  in  the  form  shown  in  Table 
7.  Ground  zero  is  at  an  altitude  of  4215  feet,  and  we  assume  that  the  sur¬ 
face  wind  was  measured  at  a  height  of  10  meters  above  the  ground.  Thus, 
the  actual  altitude  of  the  surface  wind  is  4248  feet.  Figure  9  shows  how 
cards  4a-4n  might  be  punched  for  this  set  of  data.' 


TABLE  7 

EXAMPLE  WIND  DATA  LISTING 
Site  Elevation  is  4215  Feet 


Altitude 

Direction 

(kft  relative 

Speed 

Angl e* 

_ to  MSL) _ 

(mph) 

(degrees ) 

Surface 

2 

190 

6 

15 

170 

8 

30 

180 

10 

37 

200 

30 

80 

210 

^Direction  from  which  the  wind  is  blowing. 


3 


■  . 

•i 

-•■H 

'  1 

-  H 


6.4  STORAGE  AND  COMPUTATION  TIME  REQUIREMENTS 


Approximately  15,000  (decimal)  storage  cells  are  required  by  the  CDC 
6600  computer  to  contain  the  object  codes  compiled  from  the  complete  package 
of  FORTRAN  codes  listed  in  Appendix  D. 

Computation  is  quite  fast.  A  rough  estimate  of  computation  time  is 
derived  as  follows.  In  a  single  run  on  the  CDC  6600  computer,  eight  com¬ 
plete  fallout  maps  were  calculated  and  printed,  and  for  each  map  a  vertical 
profile  of  wind  data  was  processed.  In  addition  to  the  usual  calculations, 
contour  points  were  calculated,  printed  and  punched  for  the  maps  shown  on 
pp.  58  to  80.  Activities  for  a  total  of  11,242  map  points  were  calculated, 
and  the  total  CPU  time  was  29.469  seconds.  If  we  simply  lump  all  of  the 
accessory  calculations,  including  print  and  punch  times,  in  with  the  activity 
calculation  time,  we  get  2.6  millisecond  per  activity  calculation. 


6.6  DISRUPTION  OF  CARD  INPUTS  FOR  THE  DNAf -  I  COM 


Ca  r<t 
Number 

1 


4 


4.i 


Variables  and_  format 
MOLL (12).  (12A6) 

W,  FW,  WIND,  SY,  ANG, 
GIUlFr,  (8110.0) 


ITMAP,  USETMP  (b!l) 


Cards  4a-4n  arp  input  only 
if  WIMP  -  0.0  on  card  2. 

In  that  case,  a  wind  pro- 
filp  in  (Tip  vertical  is 
specified  as  follows  (see 
sec.  6.4) 

I  OHM,  ( A4 ) 


Data  Jlescr  ip.tjon _ 

Run  identification 

W  Tota1  explosion  energy  yield  (KT) 

fW  Fission  yield  (KT) 

WIND  -  Tffective  fallout  wind  speed  (m  s  ')  (WIND  ;  0.0) 

A  blank  field  is  used  to  signal  input  of  card  4 
data . 

SY  -  Wind  shear  parameter  (s’1).  Required  only  if 
Wind  >  0.0 

ANG  Lffective  fallout  wind  vector  angle,  or  hotline 
angle,  clockwise  relative  to  north  (degrees). 
Required  only  if  WIND  ■  0.0. 

GRUFF  -  Scale  factor  for  multiplication  of  all  II  +  1  hr 

exposure  rate  values.  (See  Appendix  C.)  Default 
value  is  unity. 

II  MAP  If  true,  a  two-dimensional  H  t  1  hour  exposure  rate 
map  i s  to  he  computed  and  printed.  (See  logical 
parameter  IISLTMP  and  card  number  6.)  If  false, 

II  a  1  hour  exposure  rate  and  maximum  biological 
dose  art-  computed  for  user-specified  points.  (See 
cards  number  5.) 

MSI  IMP  If  true,  the  user  must  specify  the  fallout  map 
boundaries  and  grid  intervals  via  input,  of  card 
b.  If  false,  the  code  sets  these  parameters  via 
subroutine  SFTMP.  In  any  case,  US1TMP  is  applicable 
only  if  IFMAP  IS  TRUE. 


Specifies  whether  wind  data  are  input  in  meteorological  format 
(i.e.,  speed  and  direction)  or  in  resolved  format  (i.e.  ,  cum- 
punnnts  resolved  in  the  northerly  and  easterly  directions). 

FORM  Mi  ll  specifies  meteorological  format. 

FORM  111  SO  specifies  resolved  format 


<)l 


Card 

Number 

Variables  and  Format 

Data  Descr 

i  £t  ion _ 

4b 

FMT ,  ( 12A6) 

Wind  data  object-time  format. 

(See  catds  number  4e-4n.) 

4c 

SCALE(5),  (8F10.0) 

Wind  data  scale  factors  and  translations.  Default  values 

for  SCAFE(l)  through  SCALE (3) 

-  ).  (See  cards  4e-4n.) 

4d 

N 1 ,  N2 ,  N3,  (314) 

Wind  data  input  field  pointers 

.  Nl,  N2,  N3  are  some  permu- 

tation  of  the  integers  !,  2,  3 

(See  cards  4e-4n. ) 

4e 

AP(3),  (FMT,  see  card  4b) 

Altitude  (m  above  ground)  =  (AP(N1)  +  SCAL  1£ ( 4 )  )*SCAl.E ( 1 ) 

• 

for  FORM  MLTE : 

• 

Easterly  wind  component  (m  s  1 

) 

=  AP(N'3)*SCALL(?)*SrN(n/18D.  (AP(N2)*SCAU(3) 
t  SCALL(.)*SCALC(5>  -  180. )) 


4n  AP(N1)  999999.,  (PMT)  Northerly  wind  component  (m  s  ') 

--  AP(N3)*SCALE(2)*C0S(n/180.  (AP(N2)*SCAI.E{3) 

*•  SCAL L C 3 ) *SC At. E ( 5 )  -  180.)) 

For  FORM  RESO: 

Easterly  wind  component  (in  $"')  =  AP(N2)*SCALt(2) 

Northerly  wind  component  (m  s-1)  «  AP(N3)*SCALE(2) 

5  Cards  5d-5n  are  input  only 
if  parameter  IFMAP  ■- 
.FALSE,  on  card  3. 

5d  X,  V,  (8F10.0)  X  -  Distance  from  G7  along  the  hotline.  Positive  X  is  in 

•  •  the  downwind  direction  (in) 

Y  -  Distance  normal  to  the  hotline  (in).  Activity  and  maxi¬ 
mum  biological  dose  are  calculated  at  each  (X,Y)  point 
input. 

tin  End  of  record  flap. 

6  XMIN,  XMAX ,  YM1N,  YMAX,  Card  6  is  input  only  if  parameters  ITMAP  nid  IISETMP  both  are 

OCX,  DGY,  (8F10.0)  true  on  card  3.  Map  boundaries  and  grid  increments  are 

spec i f ied  as : 

XMIN,  XMAX  =  minimum  and  maximum  map  coordinates  relative 
to  GZ  along  the  easterly  directed  axis  (m). 

YMIN,  YMAX  minimum  and  maximum  map  coordinates  relative 
to  GZ  along  the  northerly  directed  axis  (m). 

PCX,  DGY  -  map  grid  increments  in  the  easterly  and  nor¬ 
therly  axes  directions  (in).  UGX  must  be 
specified,  but  DGY  is  optional.  If  DGY  is 
not  specified,  it  is  computed  by  the  program 
such  as  to  provide  a  spatially  und i storied  map 
when  printed  by  a  standard  line  printer. 
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To  assist  the  user  in  getting  the  code  running  on  his  own  computer, 
and  to  illustrate  the  card  input  and  the  printout  for  a  inode  1  (fallout 
map)  case,  we  provide  here  a  simple  test  problem.  The  card  input  is  given 
in  Figure  10,  and  the  complete  printout  is  given  below,  which  is  largely 
self-explanatory. 

The  raw  wind  data  are  printed  as  well  as  the  processed  data.  The 
processed  data  have  been  translated  and  scaled  according  to  the  data 
specified  on  card  4c.  Wind  layer  base  altitudes  are  in  meters  relative 
to  the  ground.  When  v,  *  and  Sy  are  supplied  directly,  this  output  is 
omitted. 

For  a  mode  2  (user  supplied  X,  Y  coordinates)  calculation,  the  output 
is  the  same  except  that  the  map  strips  are  replaced  by  a  tabulation  of 
X,  Y,  H  +  I  hour  exposure  rate  and  maximum  effective  biological  dose. 
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APPENDIX  A 


GLOSSARY  Of  SYMBOLS  AND  FORTRAN  MNEMONICS 


M(X,Y) 


FORTRAN 

Mnemonic 

OMAP(  ) , A 


DNAF1  (Function) 


CROSS  (Statement 
Fund ion) 


_ Descrl  ption _ 

li  t  1  hr  normalised  exposure  rate  (Roentgen  hr"1)  at  three 
meters  abnvp  point  X,  Y. 

Scale  factor  to  be  applied  to  activity  calculations. 

Crosswind  integrated  activity  fraction  deposited  at  distance 
X  from  ground  zero.  (m  ) 

Farfield  corrected,  crosswind  integrated  activity  fraction 
deposited  at  downwind  distance  X  from  ground  zero.  (m-1) 

Upwind  corrected,  crosswind  integrated  activity  fraction 
deposited  at  upwind  distance  X  from  ground  zero.  (m  ) 

Settling  speed  of  the  nominal  fallout  particle  at  altitude 
i .  (m  s~ 1 ) 

Factor  to  account  for  dispersion  in  the  hotline  axis  direc¬ 
tion  of  fallout  being  deposited  at  time  t  at  point  X.  (in'1; 

Activity  fraction  deposition  rote  function.  (s  ') 

Crosswind  fallout  pattern  dispersion  factor,  (nf1) 

Maximum  effective  biological  dose  (Roentgens) 

Early  cloud  radius.  (in) 

Stabilized  cloud  cap  radius,  (in) 

Wind  shear  parameter:  approximation  to  the  crosswind 
component  of  the  vertical  wind  shear,  (s'1) 


Time  of  arrival  of  the  first  fallout  at  distance  X  from 
ground  zero.  (s) 

Time  required  for  the  nominal  particle  to  settle  from  the 
stabilized  cloud  cap  base  to  the  ground.  (s) 

Yield  dependent  time  of  maximum  activity  deposition  rate 
according  to  the  q(t)  function.  (s) 

Fallout  onset  time.  (s) 
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Text 

Symbol 


FORTRAN 
Mnemon  ic 


Description 


Yield  dependent  constant  in  the  g ( t )  and  D ( X )  functions. 

(s) 

Easterly  directed  component  of  the  wind  vector  in  the 
ith  wind  stratum.  (m  s-1) 

Northerly  directed  component  of  the  wind  vector  in  the 
ith  wind  stratum.  (in  s  ') 

Effective  fallout  wind  speed.  (in  s’1) 

Effective  fallout  wind  vector.  (m  s’1) 

Energy  yield  of  the  explosion  (KT).  (The  KT  unit  is  the 
energy  released  by  explosion  of  one  kiloton  of  TNT.) 

Fission  yield  of  the  nuclear  explosion.  (KT) 

Coordinate,  relative  to  ground  zero,  in  the  easterly 
direction.  (m) 

Hotline  distance  from  ground  zero,  positive  in  the  down¬ 
wind  direction.  (in) 

Coordinate,  relative  to  ground  zero,  in  the  northerly 
direction.  (m) 

Perpendicular  distance  to  the  hotline  in  the  ground  plane. 

(m) 

Altitude  of  the  base  of  the  i'-tl  wind  stratum  hove  ground. 

(n.) 

Height  of  the  base  of  the  stabilized  cloud  cap  above  ground,  (in) 

Altitude  of  the  center  of  the  i ^  wind  stratum  above 
ground.  (in) 

Height  of  the  top  of  the  stabilized  cloud  cap  above  ground.  (ill) 

Yield  dependent  parameter  in  the  g ( t )  and  D(X)  functions, 
(dimensionl ess) 

Diameter  (uni)  of  the  nominal  fallout  particle. 

Turbulent  energy  density  dissipation  rate.  (m?  s'3) 

Variance  of  the  horizontal  spread  of  the  nuclear  cloud  at 
deposition  time.  (Does  not  include  the  wind  shear  dis¬ 
persion  contribution.)  (m? ) 

Standard  deviation  of  horizontal  spread  of  the  cloud  before 
atmospheric  transport.  (n) 
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FORTRAN 
Mnemon ic 


Description 


Wind  shear  contribution  to  the  crosswind  fallout  pattern 
variance.  (ni‘  ) 

S 1  GY  Crosswind  standard  deviation  of  the  fallout  pattern. 

ANG  Direction  angle  (clockwise  tro:u  north)  of  the  effective 

fallout  wind  vector.  (deg.) 


APPENDIX  B 


GROUND  ROUGHNESS  AND  INSTRUMENT  RESPONSE  CORRECTION  FACTORS 

To  compare  predicted  H  +  1  hour  gamma  exposure  rates  in  a  fallout  field 
with  values  measured  over  land  by  radiation  survey  meters,  it  is  necessary 
to  make  certain  adjustments  to  either  the  observed  or  predicted  values. 
Conventional  practice  is  to  adjust  the  predictions. 

Predicted  exposure  rates  are  based  on  laboratory  measurements  of  fission 
product  yields  and  on  factors  called  exposure  rate  multipliers  that  convert 
the  fission  yields  for  individual  nuclides  to  exposure  rates  at  one  meter 
height  above  an  infinite  plane  on  which  the  fission  products  are  assumed 
to  be  uniformly  distributed.  One  correction,  the  ground  roughness  factor, 
is  required  to  account  for  absorption  of  radiation  by  small  irregularities, 
or  roughness  elements,  in  an  actual  ground  surface.  The  other  correction 
is  necessary  to  account  for  variation  of  response  of  survey  meters  to  radia¬ 
tion  over  the  spectrum  of  wave  lengths  encountered.  Ground  roughness  fac¬ 
tors  for  Nevada  Test  Site  terrains  are  estimated  to  be  in  the  range  of  0.70 
to  0.75,  and  an  instrument  response  factor  of  about  0.75  is  appropriate  for 
commonly  used  survey  meters.  The  product  of  these  two  factors  is  approxi¬ 
mately  0.5,  and  this  factor  is  applied  to  dose  rates  throughout  all  of  the 
predicted  test  shot  patterns  except  for  Zuni  as  noted  in  section  5.3. 


.'i 


I. 


f” PREVIOUS  RAGE 
IS  OI.ANK 


APPENDIX  C 


FALLOUT  PATTERN  COMPARISON  BY  THE  FIGURE-OF-MERIT  METHOD 


2  0 

Rowland  and  Thompson  developed  this  method  for  comparison  of  pairs 
of  fallout  contour  maps  by  computation  of  a  single  index,  the  FM,  that  is 
a  measure  of  contour  overlap  between  them.  For  each  contour  common  to  the 
patterns,  the  area  overlapped  and  the  area  not  overlapped  is  calculated. 

The  areas  are  weighted  by  the  average  radiation  level  between  successive 
contours.  Sums  over  all  contours  of  weighted  overlapped  areas  and  weighted 
total  areas  are  computed,  and  the  FM  is  the  ratio  of  the  two  sums.  For  com¬ 
pletely  overlapped,  perfectly  matched  patterns,  FM  =  1;  for  no  overlap,  FM  =  0. 

Mathematically,  FM  is 


FM  = 


y.(ri  +  vi> 

w  2 


<ai  -  ai-l> 


A<ri  +  -V,) 

2-  2 

i  =  l 


<Ai  -  Ai-i> 


where 

N  is  the  number  of  contours  in  the  patterns.  The  summations 
are  from  highest  contour  to  lowest 

r.  is  activity  of  the  i  r.h  contour  (Roentgen  hr  *),  rQ  =  10  r^ 


a.  is  common  (i.e.,  overlapped)  area  for  the  i  th  contours. 

a  =  0. 

o 

A.  is  total  area  of  the  i  th  contour.  The  summation  in  the  de¬ 
nominator  is  computed  for  both  patterns,  and  the  largest  sum 
is  used.  A  =  0. 

o 
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The  FM  has  been  found  to  have  limited  utility  as  a  measure  of  fallout 
prediction  accuracy.  This  is  mainly  for  two  reasons.  First,  and  most 
important,  is  that  being  a  measure  of  overlap,  the  FM  is  stronqly  biased 
in  favor  of  overprediction;  that  is,  it  favors  predictions  that  cover  a 
large  area,  and  therefore  overlap  the  observed  pattern,  regardless  of 
other  considerations.  Second,  the  FM  method  imposes  no  penality  for  miss¬ 
ing  or  extra  contours;  contours  not  common  to  both  patterns  are  simply 
ignored. 
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APPENDIX  D 


FORTRAN  CODE  FOR  THE  DNAF-!  FALLOUT  MODEL 
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.  ■>  ■ 

S J8P0UT I NE  ONAF 

ONAF 

1 

n  ^ 

DNAF 

2 

c 

EXFC<JTrVE  PROGRAM  fop  THF  C-NAF-1  FALLOU"  predict  ton  model 

DNAF 

3 

c 

ONAF 

v 

c 

H,  0.  NORMENT,  AT4  0SPHE.RI3  SCIENCE  ASSOCIATES  -  AUGUST  1961 

ONAF 

5 

c 

UNAS 

6 

c 

♦  *  *  »  * 

glossary  ***»****.♦*♦.* 

DNAF 

7 

c 

A 

H+l  HOUR  :•  AMHA  EXPOSURE  RATF  (ROE  NTGENS/HR) 

ONAF 

8 

■—1  L* 

ALPHA 

onaf-i  mo:  el  paraheter 

ONAF 

9 

..  M  C 

A,T 

HEIGHTS  (•  AFOVE  NSL)  AT  HIGH  SETTLING  SPEEDS  OF  THE 

ONAF 

10 

c 

NOMINAL  PlR'ICLE  AFE  DFFINFP 

ONAF 

11 

■M  C 

AND 

HOTLINE  AJCLF  (INPUT  AS  OEGPEES  CLOCKWISF  FROM  NORTH) 

ONAF 

1? 

c 

D 

EXPONENT  ’ARAHETER  USED  FOR  APPROX.  PARTICLE  SETTLING 

DNAF 

13 

c 

SPEED  CALCULATION 

ONAF 

19 

c 

CAT 

K  FACTOR  (  (R-M** 2) / (HR-KT) ) 

ONAF 

15 

COSA 

COSINE  GF  HCTLILF  ANGLE  =  WINOY/WIND 

DNAF 

16 

€ 

OOX, OGY 

MAP  INCF.t<  ENTS  (SEE  XKAX,  ET  D  (M) 

ONAF 

17 

. m  c 

Dl  A  X 

MAXIMUM  t  r  FFCT IVF  PIOLOGTCAL  DOSE  (ROENTGENS) 

DNAF 

18 

'  .  c 

FW 

FISSION  Y [ ELD  (KT) 

DNAF 

19 

c 

GRUFF 

SCALE  FACT  OP  TO  RE  APPLIED  ~0  ACTIVITY  VALUES 

DNAF 

20 

c 

( p .  G.  GROJNC  ROUGHNESS  ♦  INSTRUMENT  RESPONSE  FACTOR) 

DNAF 

21 

c 

IFNAP 

LOGICAL  F. AG  WHICH  WHEN  TRUE  SPECIFIES  THAT  A  FALLCUT 

MAP 

DNAF 

22 

c 

31  PPEPARin 

ONAF 

23 

c  “ 

IH 

NUM8c.P,  OF  CHAPACTEF  S/IMCH  HORIIONTALLY  OH  THE  PRINTEO 

MA  P 

DNAF 

?U 

IN 

SYSTEM  IN’ UT  UNIT 

ONAF 

25 

c 

INTLIjINTl-1  FLICS  To  signal  FIRST  FASC  THROUGH  SIGTA  S  CNAF1 

ONAF 

26 

c 

10 

SYSTt  H  PR  Hr  UNIT 

DNAF 

27 

■  .  .  c 

IP 

SYSTEM  PINCH  UNIT 

ONAF 

2e 

9* 

u 

NUMBER  OF  CHAPACTEF  S/INCH  VERTICALLY  OH  THF  PRINTED 

MAP 

ONAF 

25 

c 

NHOOO 

NUMBER  OF  HEIGHTS  AT  MICH  WIND  fUTA  ARE  TN PUT 

DNAF 

3  0 

igi  R 

N1AP 

MAX.  NO.  )F  POINTS  ALLOWFO  TN  A  M A p.  (DIMENSION  CF  OMAP) 

DNAF 

31 

ll  c 

NX  MA  P 

NUMBER  OF  RAF  INCREMENTS  IN  THE  X  DIRECTION  (SEE  XMAX 

ETC) DN  AF 

3? 

c 

NY  NAP 

NUMBER  OF  MAE  INCREMENTS  IN  THE  Y  DIRECTION  (SEE  XKAX 

ETODNAF 

33 

c 

ON  A  P 

FALLOUT  Ml  P  ORDINATE  ARRAY 

ONAF 

3V 

c 

PI 

FIRL0ALL  RADIUS  ( M ) 

DNAF 

35 

c 

PS 

STARIL17E)  CLOUD  RADIUS  (M) 

DNAF 

36 

c 

SIGE 

GAUSSIAN  (AKIANCE  OF  CLOUD  riSPEPSION 

ONAF 

37 

c 

SIGY 

CROSS-WIN)  GAUSSIAN  STANDARD  DEVIATION  OF  THE  FALLOUT 

DNAF 

38 

^  c 

RATTAN  H) 

DNAF 

39 

c 

SINA 

SINE  OF  H) TL  TN6  ANGLE  =  WINDX/WIMn 

DNAF 

VO 

c 

SY 

SHEAR  PAR) METER  OR  RMS  SHEAR  PARAMETER  FROM  10  TC  GROUND 

DNAF 

VI 

c 

T  4 

TIME  OF  ARRIVAL  OF  FALLOUT  (S5 

ONAF 

42 

z 

r  ■* 

ARRIVAL  time  of  FAILCU*  FROM  TH S  STABILITED  CLOUD  CASE 

(S)ONAF 

U3 

c 

1 LL , TLS 

ARRIVAL  r 1MCS  AT  WHICH  PATE  OF  GROWTH  OF  CLOU  3  TURBULENT 

DNAF 

W 

c 

DISPERSIOJ  VARIANCE  BFCOMFS  rOHSTAHT  '") 

DNAF 

V5 

TN 

TIME  Oe  MS  XINUN  ACTIVITY  DEPOSITION  R41,;  <S) 

ONAF 

V6 

*  c 

TO 

fallout  oxsrt  time  (s> 

DNAF 

V7 

c 

USETMP 

cOGICAL  F.AG  WHITH  WHEN  TRUE  REQUIRES  THF  USER  TC  SPECIFY 

ONAF 

ve 

c 

MAP  ROUNO? IE S  AND  GRID 

DNAF 

VI- 

„  ■  ■  ■  c 

V 

SETTLING  '>  PLFDS  (M/S)  OF  THE  NOMINAL  PA  RT I  CL  E  AT  HEIGHTS 

CN  AF 

50 

c 

ALT 

ONAF 

51 

=  ’  -  c 

V  J 

settling  -peic  at  sea  level  of  the  nominal  particle  (m/si 

DNAF 

52 

c 

W 

total  yil.d  (kt) 

ONAF 

57 

HIND 

UFRCTIVr  WIND  5  PE  F  .)  (M/S) 

ONAF 

5V 

0 

W_ 

NAPFRI All  .  OGAPITHM  Or  W 

UN  AF 

55 

■  c 

ML  10 

BRIGGSIAN  10CAFITHN  OF  W 

ONAF 

56 

c 

YU 

OBSERVED  4  INC  COMPONENT'S  IN  ~HE  X  G IREC T T 0 N  (STE  /MAX 

FTC) DNAF 

57 

c 

HY 

OBSERVED  (  I  HD  COMPONENTS  IN  THE  v  DIRECTION  (S  EL  XMAX 

ETC) ONAF 

56 

c 

XNAV ,XHIN,VHAX,YM[  N  HAP  LIMITING  COORDINATES,  (X  IS  POSITIVE 

DNAF 

5° 

c 

TOWARD  LA;T  ANC  Y  IS  POSITIVE  tOWAPO  NORTH)  ( M) 

ONAF 

50 

« 


I  OH 


Pi 


.  -I 

•  i 


i 

( 

4 


■A.-.JL 


o  o  o  o  o  o 


73 
ZBH 
73  H 
70 

zr 


ST A8ILI Zc ) 

hind  laye? 

HIND  LAYER 

stabilize: 

ST  ABILI7F) 


CLOUD  BASE  HEIGHT  (  M) 
PASL  HEIGHTS  (M> 
CTHTER  HEIGHTS  (M) 

CLOUD  DFN'EP  HEIGH-  (H> 
CLOUD  TOP  HEIGHT  (Ml 


fc) 


24 


DNAF 
DNAF 
DNAF 
DNAF 
DNAF 
DNAF 
DNAF 
DNAF 
DNAF 
CNAF 
DNAF 
ON  AF 

FORMAT  (/// 15X25HG? ID  LIMITS  AND  INTERVALS/ POX 4HXNIN1 OX 4HXM AX 10X 4HYDN AF 
1MIN10X4HYMAX10X7H) ELTA  X , 3  X7HDTLTA  Y/15XF in, 0 , 4XF10. Q , 4X  Fi  0  .  0 , 4X  FI  DN  AI 


DIMFMSION  H0LM12I*  C'MAP  (  5000),  CONTUR< 
LOGICAL  I FrtAP j  USiTMP,  CFNT,  CPNC 
DATA  IN  ,  10  ,  IH  ,  IV  ,  CAY  , 

1  /  5  ,  6  ,  10  ,  6,6. 9733E9, 

DATA  IP/  7/,  S0T2J 1/2.506623275/ 


NMAP 

c000 


,  D I RC 
,57,29577951/ 


PE  ADJUSTED  TO  ACCOMOOATF  THF 


RATE 


ON  AF 
MAONAF 
DNAF 
DNAF 
DNAF 
DNAF 
(DNAF 


20,  0, 5XF10,  2.5XF10.  2) 

47  FORMAT <//15X,  61H)GX  AND  3CY  HUS 
IP  IN  STORAGE) 

1000  FORMAT  (8F10.0) 

1100  FORMAT ( 5L1) 

12  00  FORMAT  ( -*  (  5X,  IPEll.u)) 

1300  FORMAT  (///7X,  1181A  MAP  IS  NOT  PREPARED.  H*1  HOUR  EXPOSURE 
1R0ENTGENS/HR)  AND  MAXIMUM  EFFECTIVE  BIOLOGICAL  DOSE  (NEBO) 

2EN  S ) /  tx,  J6HARF  lOMFUTEO  FOR  USER-SPECIF IFD  POINTS.  X  IS 
30WAP.D  DIRECTION  AMD  V  IS  OF  OSSWIND,  // 

4  9  X ,  hHX(M),  i?X,  4HY (M)  ,  11X,  6HR(H+1),  11*,  4HHEB0) 

14  00  FORMAT ( 1 2 A6) 

15  00  FORMAT  (  1H1,  50  <  ,  19H*  ,**.**♦*  *//55X,UHD  N  A  F 

1  2  5  X ,  7  2HT  H  E  )  E  E  E  N  S  E  NUCLEAR  AGENCY 

2T  F  A  L  L  0  U  T/'  4§X,  3?HP  P.EDICTION 

3  E  M//  51X,  ll?K*  M  *•»*»*.  *■//// 65  X ,  11H  PREPARED  BY/46X,  30HQNAF 

4  AT  MO  SPH  ERIC  SCIENCE  ASSOCIATES/  54X,  14HBFDFORD,  NASS . ////2 5X ,  DNAF 
5 30  H****  RUN  IOENr IFICAT I3N  **•*,  3X»  12A6) 

16  00  FORMAT  (//20X,24HYI  CLDS  -  TOTAL  .FISSION),  21X,  1PE12.5,  2H  (, 

1  IP F 12. 5,  4H)  XT/  20  X  , 21  HO ROINATK  SCALE  FACTOR,  ’4X,  1PE12.5) 

1700  FORMAT (//’OX,  20HIFFFCTI VE  NINO  SPFED,  25X,  1PF12.5,6H  M/SEC) 

18  00  FORMAT  (  ’OX,  37HHI  t  ID  DIRECTION  (CLOCKWISE  FROM  NORTH),  8X, 

1  IPE12. 5,  5H  PEG,  i 

1  9  0  0  FORMAT ( 20X,  19HSHCAR  PAR  AS  F TER,  SY,  26X,  1PE12.5,  an 
20  00  F.'RMAT  (//20X,  71H5M  EFFECT  I VE  WIND  SPEED  LESS  THAN  0 
1EN  SPECIFIED  -  WI«ID  =  1PE12.5/  20X,  RxHWINr  IS  RESET 
2100  FORMAT ( Ab) 


(ROENTGQNAF 
THF  WINDNAF 
ON  AF 
DNAF 
DNAF 
DNAF 
DNAF 
TDNAF 


-  1// 
FAS 
SYS 


PER  sro 
.50  M/S  HAS 
TO  0.50) 


CR0SSWIN3  G) US  SI AN  SPREAD  VIA  THIS  FUNCTION 
*Y,SIGYY)  =  FXF (-0.5*  ( vy/S IGYY> 2) /$ IG YY/SQT2P 1 


HO  1 L 


COMPUTE  THE 
CROS 

c 

COPY  IN  RUN  IDENTIFICATt ON 
100  READ  (  IN,  l-*00,tND=;  00) 

C 

150  INTL 1=0 
1ST  L  2=  0 

COPY  IN  BASIC  DATA 

R- AD  (IN,  10011)  W,  FW, 

READ  (IN,  llOO)  IF-IAP, 
wlio  =  alogio(w) 

ANO  -  AN5/DIPC 
SINA  =  3  IN ( AND) 

COS  A  =  COS(ANG) 

IF  (GRUFF  .  EQ.  O.m  GRUFF  =  1,  D 
C.AYFW  =  CAY*FW‘GRJFF 

CO-FiJTF  STABILIZED  CLOU?  PROPERTIES 


WIND,  SY, 
US  El M  P 


ANCi,  C.PUFF 


DNAF 
DNAF 
DNAF 
DNAF 
DNAF 
ON  AF 
DNAF 
BEONAF 
DNAF 
ON  AF 
DNAF 
ON  AF 
DNAF 
nN  af 
ON  AF 
DNAF 
DNAF 
DNAF 
DNAF 
DNAF 
DNAF 
ON  AF 
DNAF 
DNAF 
DNAF 
CNAF 
DNAF 
DNAF 
DNAF 


61 

6? 

63 

64 

65 

66 
67 
6E 

69 

70 

71 

72 

73 

74 

75 

76 

77 
76 
70 
80 
81 
8? 
63 

84 

85 

86 
87 
86 

89 

90 

91 
9? 

93 

94 

95 
9b 

97 

98 

99 
100 
101 
102 
1  03 

104 

105 
10b 
107 
lOfc 

109 

110 
111 
11? 

113 

114 

115 
life 
117 

116 
IIP 
120 


109 


■A 


*  , 


7S 

■ 

a 

a 

-d 

Ml. 


CALL  SCLOP(  W,  RI,  RS,  ZT,  ZB) 

DNAF 

121 

COPY  0 LIT  BASIC  DATA 

ON  AF 

12? 

WRITE(IO,  1500)  H3LL 

DNAF 

123 

9 

WRITEdO,  1600)  Hi  FW,  GP.JFF 

DNAF 

124 

IFCHIND  . Me.  0. 0)  GO  TO  ZOO 

DNAF 

1  25 

COPY  IN  WINO  PROFILE  DATA 

DNAF 

126 

CALL  INNINOIIN,  13) 

DNAF 

127 

COO  PUT  £  EFFECTIVE  WIND,  HIM3 

DNAF 

126 

CALL  EFHINOf  ZT ,  ZB,  WIND,  SIIIA,  COSA) 

ON  AF 

125 

COM PUTf  WINO  SHEAR  PARAMETER,  SY 

ON  AF 

130 

r 

SY  •-  SYHNO  (ZT,  T 3,  SINA,  COSA) 

DNAF 

131 

CHICK  WINO  FOR  THRtSHOL;  AND  SET  WINr=0.50  IF  1“  IS  BELOW 

THE 

THRESHOLD 

DNAF 

13? 

200  IF  (HIND  .Grl.  0.50)  GO  TC  20  5 

DNAF 

133 

WRITEdO, 2000)  WHO 

DNAF 

134 

WINO  =  0.50 

DNAF 

135 

COMPUTE  FALLO'JT  ONSET  TIME,  TO 

DNAF 

136 

205  TO  =  ONSET (W) 

DNAF 

137 

COPY  OUT  MODEL  P4RAMFTEJS 

DNAF 

13f 

9 : 

WRITEdO,  1700)  WING 

DNAF 

135 

j 

IF  (  .NOT,  IFMAP)  !-0  TO  210 

DNAF 

140 

» 

DIR  =  OIRC  *  ATAN< SIMA/COSA  ) 

DNAF 

141 

IF  <  COSA  ,LT.  0.0)  DIR  =  DIP  -  SIGNdflO.O,  IIP) 

DNAF 

14? 

WR ITS ( 6 ,  IB  00)  DI? 

DNAF 

143 

210  WRIT£(IO,  1900)  Si 

DNAF 

144 

• 

CHECK  IF  A  MAP  IS  TO  Bt  PREPARED 

DNAF 

145 

r-; 

IF ( I FMAD )  GO  TO  300 

DNAF 

146 

HR  IT  E ( 10  ,  1300) 

ONAF 

147 

■_  *1 

COPY  IN  COORDINATES  AT  (HICH  ACTIVITY  is  "0  BE  COMPUTED 

DNAF 

1  46 

*  ■  , 

250  RE AO ( IN, 1 OOfl ,EMD“1  00)  Xj  Y 

DNAF 

149 

W  ^ 

C 

DNAF 

150 

COMPUTE  CLOUD  DISPERSION  PARAMETERS  AND  FALLOUT  ARRIVAL  TIME 

ON  AF 

151 

275  CALL  S1GT  A<  <  ,  W,  TO,  HIND,SY,  P.I,  RS ,  ZB, 

7T, 

WL  10, 

DNAF 

152 

it! 

1  TNT1.1,  SIGH,  SHY,  TA) 

DNAF 

153 

COMPUTE  CROSSWINO-INTEGLATtD  ACTIVITY  FRACTION  ON  THE 

HOTLINE 

AT  X 

DNAF 

154 

GO FX  -  ONAFKX,  MIND,  SIGE,  H,  INTL?) 

DNAF 

155 

COMPUTE  H  +  l  ACTIVITY  <  R'  HR)  Ar  POINT  (X,Y) 

DNAF 

156 

A  =  CAYFH  -  CROSSC Y.SIGY)  *  GOFX 

DNAF 

157 

COMPUTF.  MAXIMUM  EFFECTIVE  IIOLOGICAL  OOSE 

ONAF 

156 

OMAX  =  A  *  OOSF  <  T  S ) 

DNAF 

155 

'ml 

COPY  OUT  RESULTS 

DNAF 

160 

WRITEU9,  1200)  X,  Y,  A,  .  ,-!AX 

ONAF 

1  61 

CO  TO  250 

DNAF 

162 

r..  • 

C 

DNAF 

163 

r.'  - 

CHECK  IF  MAP  IS  SPECIFIED  GY  rHE  USFf  OP  Bv  THF  PROGRAM 

DNAF 

164 

.V 

300  IF(USETMR)  GO  TO  *00 

DNAF 

165 

CALL  ScTMP(W,  FH,  WIND,  SINA,  COSA,  GRUFF, 

XMIN, 

XMAX, 

ONAF 

lb6 

'■  • 

1  Y  MIN,  YMAX,  DGX,  DGY) 

DNAF 

167 

St, 

GO  TO  *,50 

DNAF 

166 

A 00  READ  (IN,  10  00)  XMI  N,  XRP  X,  V*Ii:,  YMAX,  OGX,  DGY 

ONAF 

169 

// 

410  IF  ( DGY  0,0)  3  GY  =  DGX* IH/T V/2, 0 

DNAF 

170 

450  NXHAP  ~  (XMAX  -  K1 IN ) /DG  X 

DNAF 

171 

NYMAP  a  (YMAX  -  V«IN)/r>CY 

ONAF 

17? 

MM AP  =  NXMAP*NYNA3 

ONAF 

173 

IF  f  MMAP  iLE-.  NMAM  GO  TO  460 

ONAF 

174 

WRITE(I0,47) 

DNAF 

175 

DGX  =  SQRT(  ?.  0* 1/ * MMAP* DGX* DGY/ <NMAP*IH) ) 

ONAF 

176 

OGY=0. 0 

DNAF 

177 

‘ '  - 

CD  TO  -»10 

DNAF 

176 

4  bQ  WRITE  (10,  2h)  ;<Mt  N,  XMAX,  THIN,  YMAX,  OGX,  DGY 

ONAF 

1  7° 

COMPUTE  ACTIVITY  AT  MAP  POINTS 

ONAF 

160 

Jb 


L 


1 10 


VI  =  VMIN 

ON  AF 

1S1 

INDFX  =  0 

ONAF 

187 

DO  500  J=1,NvmAP 

DN  AF 

ie3 

YM  =  YM  +  OGY 

ONAF 

1  84 

XM  -  XMIN 

ONAF 

185 

DO  500  I=1,NXM«F 

ONAF 

186 

XI  =  XM  >  OGX 

ON  AF 

187 

IF  (NINO  .GT.  l.OE-02)  C-D  TC  *t70 

ONAF 

186 

X  =  SORT CXM**?  *  f  M»* ?) 

ONAF 

189 

Y  =  0.  0 

ONAF 

190 

GO  TO  47? 

ON  AF 

191 

470  X  =  XH+SINA  >  YD^IOSA 

ONAF 

19? 

V  s-XM^COSA  +  Y(:*»  ICA 

ONAF 

193 

CCMPUTF  CLOUD  OISRFPSIGV  PARAMETERS  AND  FALLOU 

T  ARRIVAL  TIME 

ON  AF 

194 

>72  CALL  SIGTAdr  H »  TO,  WIND, SY ,  PI, 

P.S,  70,  7T ,  WLIO, 

ONAF 

195 

1  INTL1,  SIGE,  SliY,  TP) 

DNAF 

1 9f 

INDEX  =  INDEX  >  l 

DN  AF 

197 

COMPUTE  MAP  ORDINATE  -  M  +  l  HDUP  EXPOSUPF  RATE 

<  NORMAL  1 7ED) 

ONAF 

196 

475  OMAP ( INDEX)  -CA yf* * CNAF1 < X ,  WIND,  STGE, 

W,  INTL2)  *  CROSS! V , SIGY ) 

ONAF 

1  99 

500  CONTINUE 

ONAF 

200 

CO?Y  OUT  MAP 

ONAF 

?D1 

6  00  CALL  M  APP  (  ID,  HOLL,  DGX,  nr,Y, 

NXt'AP,  NWMAP ,  OMAP, 

ONAF 

■>o? 

1  XHAX ,  XMIN,  Y  H  A { ,  VM1F) 

DNAF 

203 

GO  TO  100 

DNAF 

Z04 

900  RETURN 

ONAF 

?06 

END 

DNAF 

2  Db 

SUBROUTINE  SCLOf' <  W,  Rl,  RS,  7~,  IB) 

0 

COMPUTES  STfteiLIZSU  CLOUT'  PROPERTIES 
C 

C  H,  G.  NORM  ENT ,  ATMOSPHERIC  SCIENCE  ASSOCIATES  -  AUGUST  19ED 

C 

RI  =  1  91.  <9.  S3) 

WL=ALOG ( W> 

RS  -  EX°(o.7553  +  WLMO  3?05r>  +  WL*  0.  0  i  1  3?4  78 )  ) 

IF  (M  .GT.  >.07)  CO  Tf  10 
ZB  "  ?22  i.  *  W**  ( 01  31*63) 

GO  TO  10 

10  7R  =  166 1.  *  WMO.BlSe) 

20  IF  <  H  .GT,  ?.29>  IF ( W  -  19.  0)30,30,40 
7T  -  35)7.  *  H**  ( 0*.  ?5  53) 

CO  TO  60 

30  ZT  =  3170.  *  H**(n.L077) 

no  to  so 

40  ZT  =  G>74,  »  W^d'.lf  50) 
so  continue 
RETURN 
END 


SCLOD  2 
SCLOC  3 

scloc  4 

SCLOC  5 
SCLOC  6 
SCLOD  7 
SCLOr  6 
SCLCO  9 

scloc  in 

SCLOD  11 
SCLOC  1? 
SCLOD  13 
SCLOC  14 
SCLOO  15 
SCLOO  IF 
SCLOD  17 
SCLOD  16 
SCLOO  19 
scLon  2n 
SCLOO  21 
SCLOD  2? 
SCLOT  ?3 


OOOOOOOOOO 


SUBROUTINE  INWIFDI  IN,  10) 

0BI0S  IN  AND  PROCESSES  OBSERVE:)  WlhJf  OAT  A 


H,  G.  NORHENT,  AHOEPHLf.IC  SCIENCE  ASSOCIATES  -  AUGUST  I960 

4  44¥¥4.¥-444444444444-.4  »¥**.**4*44  44*444 

RE  60S  ANJ  PROCESS:  5  HIND  3/TA  FOR  A  HORIZONTALLY 

HOMOGEN IOUS  FIELD.  VER'ICAL  COMPONENTS  ARE  NOT  CONSIDERED. 

C  **4.*»»»*.*t|»»  GLOSSARY  *4.».-4»4¥4>4444» 

C  FMT  0<) JEC’-TI-'E  FORMA!  OF  WIND  DA~A 

C  FORM  INPUT  PARAMETER  to  INDICATE  ECRHA"  OF  HIND  DATA  TO  FOLLOW 

C  •EITHER.  -  3FTE0R  OR  RESOL V 

C  N1,N2,N3  DATA  El: LD  POINTERS 

0  SCALE  DATA  SCALE  r ACTCRS  AND  TRANSLATIONS 

C  SEE  GLOSSARY  IN  P!  OCCAM  O')  A  E  FOR  OTHER  OUANTIES 

g  4  444  4444  444  4  44444444  4*  4  44 444 4*4 44**444444444«« 4444444.41 4* 44 *1  •  444444 

REAL  METEOR 

COMMON  /WHAT/  NH0)0,  7PH(50),  7CH<cQ),  WX(FD),  WY(50> 

DIMENSION  SCALE  <  )),  AP(  3),  FHT<1?> 

C 

OATA  ALIMIT  ,  RA C  .  PROGRM  ,  HrTEOF  ,R£SPLV 

1  f  999999.  ,  01r  4532925,  6HTNHIND,  4HHETE  >  4HRFS0/ 

OATA  IREC/62 

C 

l  FORMAT (  4X,  6HLE ViLS,l4,  5H  "HRU,I4,  0F12.5) 

3  FORMAT  (  7///33X,  »  5HWIND  LAYER  BASE  ALT  fUCES/) 

C  ¥  FORMAT  <  1H03X31HMA<  IF.UM  WIND  SPACF  ALTITUDE  IS  E17.5, 7H  METERS) 
1000  FORMAT  (T2A6) 

UOU  FORMAT  (3F10.0) 

1200  FORMAT  (2014) 

13  00  FORMAT (////  26X,  13HRAH  HI  MO  DATA, TTX, 19HER0CESSE0  MING  DATA/18X 

11HZ,  9X,  10HVX  CR  DTP.,  -3X ,  11HW  OP.  SPEED,  14X,  1H7,  12X, 

2  2HW  X,  1»X,  ZHVY) 

1400  FORMAT  (10X,  3(2X,  1FE12.5)) 

15  00  FORMAT  <1H+,57X,  M2X,  1PE12.5)) 

IS  00  FORMAT  (  1H 0 ,  EX,  I4HWIMT  STRATA  ALTITUDES  INC ONS TSTE NT ) 

1300  FORMAT(Ah) 

0 

CJ3Y  IN  DATA  SPECIFICATION 
READ  (IN, 1300)  FOR* 

CHECK  FORM 

IF  (FORM  .EO,  meteor.)  c-C  TO  75 
20  IF  (  E  ORM  .ME.  PESO. V) CALL  E P ROR < PROGRM, -70 , TO) 

C03Y  IN  FORHAT,  SCALF  T,  FITLD  POINTERS 
25  P.EAO  (  IN  ,  10U0)rMT 
READ  (  IN  ,  LlOui  SCALE 
READ  (  IN  ,  1200)  N1 ,  N2  ,  (13 

no  so  i  -  1,3 

50  IF(SCALLU)  .EQ.  ’  0.  0)  SCALE  <  I  >  =  1.0 

IF  (  FORM  .EQ.  MFTF: R)  T PUS* S C ALE ( 5) * SCALE ( 3)  -  lg 0. 

WRIT-  (10, non) 

MHO!)  0  =  0 

C  0 3  Y  IN,  PRINT  RAW  DATA,  TRANSLATE  AND  SCALE  DATA,  ADD  PRINT  PROCESSED 
C  DATA 

100  READ  (  IN  ,  KMT)  VF 

IF  (AP(Nt).  GE,  ALIMI  T)  CO  T(i  ?5  o 


INWIN  2 
INWIN  3 
INWIN  4 
INWIN  5 
INHIN  6 
INWIN  7 
• IN  WIN  6 
INHIN  9 
INHIN  10 
INHIN  11 
INHIN  12 
INWIN  13 
INWIN  14 
-INWIN  15 
INHIN  It 
INWIN  17 
INHIN  16 
INHIN  IP 
•INHIN  20 
INHIN  21 
INHIN  22 
INHIN  23 
INHIN  24 
INHIN  25 
INHIN  26 
INHIN  27 
INHIN  26 
INHIN  29 
INHIN  3C 
INHIN  31 
INHIN  32 
INHIN  33 
INHIN  34 
INWIN  35 
, INHIN  36 
INHIN  37 
INHIN  36 
INHIN  39 
INHIN  40 
INWIN  41 
INHIN  42 
INHIN  43 
INHIN  44 
INHIN  45 
INWIN  46 
INHIN  47 
INHIN  46 
INWIN  49 
INWIN  50 
INHIN  51 
INHIN  52 
INHIN  53 
INHIN  54 
INHIN  55 
INHIN  56 
INWIN  57 
INHIN  56 
INWIi  59 
INHIN  60 
INHIN  6 ) 


WRITE!  10,1400)  ftPHl)  ,  AP(N2),  AP<N3) 
mODO=NHOOO«-l 

7CH<NHOOQ)=!AP(N1I  +  SCC L£ ( l> ) *SCALE!1) 

IF ! FORM. £Q. RESOLD  GO  TO  150 

NX<NHODO>  =  APOJ3>  *SCALF  1 21  *  SIN  I R ACC* ( A  PI N2 > *SCAL £ !  3)  +  TRNSJ) 

HY(NHOOO)  =AP!N3)»SC6LF(2)JtCOStRADC,*<AF(N2)*SCALe<3)  *  T  RMS) ) 
GO  TO  200 

150  WK(NHOOO)  -  >  P!N2)*SCAC£<2» 

WYCNHODQ)  =  5P(N3)'SCAL£!?) 

200  WRITEIIO,  1500)  ZC! < NPODO ) , W X INHODO) ,  WVfNPODO) 

GO  TO  100 

COMMINGLE  DATA  TO  ARYANS  E  IT  IN  CROEP  OF  ASCENDING  ALTITUDE 
250  NHODMl=NHOOO-l 

DO  255  I=1,NHOOM1 
IP1=I*1 

DO  255  J-IP1,NH0D! 

IF!  ZCH(I)  .LE.  ?  CH 1 J) )  GO  Tfl  z 55 
TE  NP=  2CH  (I) 

2  CH ( I ) =  ZCH  < J) 

ZCH  < J) =T£MP 
TEMP=HX(I) 

WX  <  I )  =WX  <  J) 

WX  < J) =  TENP 
T£MP=WY!I) 

WY  <I)  =  HY<J) 

WY(J)=TEMP 
255  CONTINUE 

CONSTRUCT  HIND  LAYER  65 SE  ALTITUDES  IN  ARRAY  ZBH 

259  Z3H  < 1)  =  0.0 

00  260  I=2,NH000 

260  ZBH(I)  =  (ZCH(I-l)  *  ?Ch(I))/2,  9 

C  ZNAX=2. Q+ZCH(NHO0))  -  ZC H(NHODO) 

COPY  OUT  WIND  LAYER  8AS:.  CATA 

WRITE!  10  ,3) 

DO  270  IG0=ljNH001 > IREC 
ISTOP=IGO+IREC-l 
IF  1 1  STOP  .GT. NHOPOI  IST0P=N40D0 
270  WRITE!  10  ,  1)  IC-O,  ISTOF,  <  ZP  H(K)  ,  K=IGO,  ISTOP) 

C  WRITE!  10  ,A>ZHA< 

RETURN 

FND 


INWIN  6? 
INWIN  63 
INWIN  EL 
INWIN  6E 
INHIN  66 
INWIN  67 
INWIN  66 
INWIN  69 
INWIN  70 
INHIN  71 
INWIN  72 
INWIN  73 
INWIN  7L 
INWIN  75 
INWIN  76 
IN  KIN  77 
INWIN  76 
INWIN  79 
INWIN  eo 
INWIN  81 
INWIN  82 
INWIN  83 
INHIN  64 
INWIN  85 
INKIN  86 
INWIN  87 
INWIN  86 
INWIN  89 
INWIN  90 
INWIN  91 
INHIN  9? 
INHIN  93 
INWIN  94 
INWIN  95 
INWIN  96 
INWIN  97 
INWIN  96 
INWIN  99 
INWINIOO 
INWIN101 
IN  WIN  102 


SUBROUTINE  FTWINPdT,  7B,  HIGH,  SINA,  COCA) 


EFWIN  2 

C  FFWIN  3 

COMPUTES  THE  ;FFICTrV£  W I  hip  SPEED,  WIND,  AND  tTS  PIRZC’ICH  ANGLE  FROM  EFWIN  4 

C  OBSERVED  HIND  nA  Ti1  E  F HI  N  5 

C  FFWIN  6 

C  H,  G .  NORMENT,  ATMOEPWEPIC  SCIENCE  ASSOCIATES  -  JULY  1931  EFWIN  7 

C  EFWIN  o 


COMMON  /H.JAT/  NHOOO,  7PH<5(I),  ’CH(5n>,  W<<F0>,  WV(5D) 

DIMENSION  ALT(  39)-,  V(  39) 

DATA  ALT/0  ,  ,  100  0,  ,.20  00t  ,  3003.',  <-100.  ,  5100,,  6000.,  7  000,,  «000 


1 

5ono. , 

10000., 

1  1)0  0  0., 

1  20  00.  , 

13000.  , 

14000.,  \ 

5000.  , 

16000.  , 

7 

17000. , 

1  3000.  , 

isnoOl , 

2  00  10  , 

21000  , 

’2000. , 

23000., 

24000.  , 

3 

25000., 

26  0  00.  , 

7?  1)0  0.  , 

30000., 

72  0  0  0.  , 

34000., 

36000., 

3e  aoo., 

*♦ 

40  0  0  0., 

•12000.  , 

-4000,  , 

4  60  00.  , 

42000  , 

50000./ 

DATA  V/ 

1.6636, 

i.n?4. 

1.  11  7h4  > 

1.6401 , 

1 «  90  97,  ?. 

.96  36  , 

2.  0621, 

1 

2. 1451, 

7. 2.35 n. 

7.  3303, 

2, 4  3?  4, 

2  5419, 

2,  6446, 

2.7469, 

2.  6  551  , 

2 

?. 9630, 

3. 0722, 

ME  31, 

3.2951, 

7.403’, 

3  .5222, 

3.6322, 

3.  P 122, 

3 

3 . 9  2*,  7  , 

7.999-, 

oeoc, 

4. 1 ’46, 

4.  76  65 , 

4,  5326, 

4.7063, 

4.  651  7, 

4 

6.  0  023, 

6.1767, 

'.3910, 

5.6672, 

6.0155, 

6.4727, 

7 . 0  7  7  6  , 

7.  6819/ 

c 

COMPUTE  CLOUD  CENTER  HL  If-PT ,  70 
70  =  { ZT  +  70)/?. o 
COMPUTE  LOCATION  OF  70  UN  ZBH  A  RPAV 
MHODO  =  KiHOOO  +  1 
no  100  I=l,NHOno 
J  =  MHODO  -  I 

IF  (  Z  0.  GT  ,  78H  (  J )  )  GC'  TO  15  0 
100  CONTINUE 

150  7  =  <ZC  +  Z3H(J))/2.f 

CALL  T R 3 L ( 7 ,  35,  °LT ,  V,  VP) 

WGT  =  <Z0  -  79H<J))/VP 
SWOT  =  HGT 

CALL  T  R  0  L  <  7 ,  NHOOJ,  7CH ,  W  X, WNTX) 
CALL  TRdL(2,  NHOCJ,  7CH,  WV,WNPY) 
WNOX  =  WNOX* WGT 
HNOY  =  WN)Y*WC,T 
J  =  J  -  1 

I F  <  J  .  EQ.  0)  GO  TO  2  5  0 

no  zoo  i=i, j 

CALL  T3’(.(  7CH  ( I )  ,  36,  ALT,  V,  VP) 
HGT  = ( 2  3  H  ( I *1)  ~  7DH ( I) ) /V  P 
swgt  =  shgt  ♦  hot 


EFWIN  c 
FFWIN  in 
..FFWIN  11 
EFWIN  1’ 
EFWIN  13 
FFWIN  lu 
EFWIN  IF 
EFWIN  16 
EFWIN  17 
FFWIN  IP 
EFWIN  1<= 
EFWIN  ?0 
EFWIN  ?1 
EFWIN  27 
cFWIN  23 
FFWIN  ?t 
EFWIN  25 
EFWIN  ?b 
EFWIN  27 
FFWIN  26 
EFWIN  2C 
EFWIN  30 
EFWIN  31 
EFWIN  3? 
EFWIN  33 
FFWIN  3 k 
EFWIN  36 
EFWIN  36 
EFWIN  37 
EFWIN  36 
EFWIN  39 
cFWIN  4(1 
EFWIN  61 
EFWIN  6? 
FFWIN  U3 


HNOX  = 

WNOX  4-  WX(I)»WCT 

EFWIN 

44 

20  0 

WNOY  = 

WNOY  +  WY  Cl ) *  WGT 

EFWIN 

45 

250 

WIND  = 

SORT (WNnX**?  ♦  WHO y*M2) 

EFWIN 

4f 

IF ( WIND 

.FQ.  0.0)  C-C  TO  3  0  0 

FFWIN 

47 

SINA  - 

WNOX/WIMD 

EFWIN 

46 

COSA  - 

WNOY/WIND 

t  FWI  N 

4C 

WIND  = 

W  INO/SWG"r 

EFWIN 

5  0 

RETURN 

FFWIN 

51 

3on 

SINA  s 

1.0 

EFWI  N 

5? 

COSA  = 

0.1 

EFWIN 

53 

RETURN 

EFWIN 

54 

END 

FFWIN 

5r 

FUNCTION  DOSE.  (TA ) 

c 

COMPUTES  TH£  FACTOR  WHICH  WHEN  MULTIPLIED  PY  A  PROVIDES  MAXIMUM 
C  EFFECTIVE  eiOLOCIItL  DOTE 

IF  ( T  A  .LT.  UP7.9i  GO  TO  100 

71E  --  2.0  *  ALOGtr A/3600. 0)  +  4,0 

OOSE  =  4,618?  -  2:'F  *  (0.53567  -  0.  0169?3  *  ZEE) 

RETURN 

100  ZEE  =  Q.8685S33  *  A  L  C'  G  <  T  A ) 

OOSE  =  15.?891  -  TEL  *  (?.903??5  -  0.166?3?15  *  ZEE) 

RETURN 

ENO 


DOSE  1 
DOSE  ? 
DOSE  3 
DOSE  4 
DOSE  5 
DOSE  6 
DOSE  7 
OOSE  6 
DOSE  5 
DOSE  10 
DOSE  11 
DOSE  1? 
OOSE  13 


FUNCTION  SYWNO(ZT„  ZD,  FIMA,  COSA)  SYHND  ? 

4  SYHND  3 

COMPUTES  RMS  SHEAR  PARAMETER  PE  T  KEEN  ZO  AND  tH-E  FROUND  SYHND  4 

C  SYHND  5 

C  H.  G.  N3RMENT ,  ATMOf PHEPIC  SCTFNCE  ASSOCIATES  -  AUGUST  1S80  SYHND  6 

c  SYHND  7 

COMMON  /  HTAT/  NHODO,  ZDH(SO),  "'CHI  5  0),  WX(50),  HV(50I  SYHND  f 

COMPUTE  CLOUD  CiNT_R  HEI'GH-  AI  D  HALE  CLOUD  THICKNESS  SYWNC  9 

70  =  <ZT  F  79)/?.  O’  SYHND  10 

DM  =  (ZT  -  ZD)/?,  0  SYHNU  11 

N  =  0  SYHND  1? 

S*  r  0.0  SYHND  13 

7U  =  70  ♦  OH  SYHND  14 

7L  -  ZO  -  OH  SYHND  15 

COMPUTE  SQUARE  OF  SHE AP  RAF  AMLTE  P  AND  INCREMENT  SUM  AND  COUNTER  SYHND  16 


CALL 

TRRL(7U,  NHODO,  zee, 

HX , 

HXU) 

SYHND 

17 

CALL 

TRPL(TU,  NHC'ZO,  7CP, 

HY, 

WYU) 

SYHND 

18 

CALL 

TRPH7L,  HHOOO,  ZCH, 

HX, 

HXL) 

SYMNO 

19 

CALL 

TRRH7L,  MHOCO,  7CH, 

HY, 

WYL) 

SYHND 

20 

SY  = 

SY  ♦  ( ( -COSA *  (WXU  - 

H  XL ) 

+  SI N  A» ( W YU  -  WYL)  )  /  (ZU-  ZD)  *  •? 

SYHND 

?1 

N  =  N 

♦  1 

SYHND 

2? 

IF  DONE 

SYHND 

?3 

IF  (ZL 

•EO,  0.0)  G j  TO  ?00 

SYHND 

?  4 

7U  = 

7L 

SYHND 

?5 

7L  = 

7  U  -  ■».■)*•  OH 

SYHND 

26 

CHFCK  IF  THIS  IS  LAST  SUMf'A  NT 

I F ( 7  L  .LT,  7DH(?))  7L=0.0 

GO  TO  100 

?00  SYHND  =  SDRT(SY/N) 

°  F  T  U  R  N 
F  N  0 


SYHND  ?7 
SYHND  28 
SYHND  ?9 
SYHND  30 
SYHND  31 
SYHND  3? 


o  o  o  o  c>  o 


FUNCTION  'INS -IT  (W) 

C 

COMPUTES  FALLOUT  ONSiT  T'IPL  (?) 

C 

C  H.  G.  NORM -NT,  ATMOSPHERIC  SCIENCE  ASSOCIATES  -  AUGUST  198T 

L 

if(h  ,lt.  i.o.io  so  To  100 

ONSFT  =  1147.54 
RETURN 

100  WL  -■  ALOG(W) 

ONSrT  -  FXP(  l.c275<S7  +  WL»  (  0.  403<5<.66  +  HL*  ?«  0  643? 2E- 2)  ) 

RETURN 

END 


ONSET  2 
ONSET  3 
ONSET  4 
ONSET  5 
ONSET  6 
ONSET  7 
ONSET  8 
ONSET  g 
ONSET  10 
ONSET  11 
ONSET  1? 
ONSET  13 
ONSET  14 


SUBROUTINE  SIGTA(<,  W,  TO,  WINO.SY,  PI,  PS,  7B,  7T,  WLIO, 

1  INTL1,  sige,  shy,  tm 

OiPUTES  GAUSSIAN  VARIANCE  SIGE  AND  STANDARD  DEVIATION  Si'GY,  AND  TIME 
OF  ARRIVAL  OF  FAL.OUT  AT  DISTANCE  X  FROM  GZ  ALONG  THE  HOTLINE 

H.  G.  NORMcNT »  AT *  OSFHtR IC  SCIENCE  ASSOCIATES  -  AUGUST  1501 

DATA  GX  VO  ,  R  ,  EX 

1  /  *♦.  7D607E-5  ,  2.9E-5  ,  0.66666566667  / 

C 

IF ( INTL1  .GT .  0)  10  TO  100 

q  ♦  »»****♦******#•<»•*»*»».********«  »»* 

COMPUTE  INITIALIZATION  faf  amltep.s  on  first  PASS 
T3  -  (1.0  -  EXP (-3 *ZR) J/BXVO 
wro  =  wino  *  to 
TO PP.  IW  =  TO  +  P.I/tlhP 
RSO?  =  RS/7 . 0 

IF ( ML  1 0  ,GT.  1.0)  IF(HL10-3.0) lfl, 20,20 

SI  GO  P  =  RI 
GO  TO  ?Q 

10  SI  GOP  =  RI*  <  1.  0  f  3 1  P*WL  10) /«.  0 
GO  TO  SO 

?0  SIGOP  =  2.5*RI 

50  CINT  =  (WSOZ-SIGO5  )/<TP  -  ’0) 

SIGC1S  =  SIGOP**E( 

SIGC1L  -  RS02«MX 

SIGC2S  =  T.oF.e*si;  ci?  -  ?.  oeb 

SIGC2L  =  T.0E6»Si;  C1L  -  2.0E9 
E=SC1  -  0-  0 165  ?  ?*  ■<  ••  <  -  0.  10  2  33) 

EPSC2  =  T.  QF.Cr'FTS:  1 

t.s  =  (imo.n  -  st  gcis)/fpsci 

TLL  =  (1000, n  -  St GC1L)/£DSC1 
SH2ARC  =  SY*(ZT  -  ZEO/IO.O 
INTL1  =  1 

.*»»♦»».»»****♦*■»****♦******* 

COT  PUT  E  FALLOUT  TIME  OF  ARRIVAL 
100  XI  Ml  =  X/WINO  -  TIPFIH 

TA  -  TO  +  (XIM1  F  SRFKX  IMl-»*2  ♦  0 . 0  l*TOPP  IW)  )  /  ?,.0 
COMPUTE  CLOUD  HORIZOnTA..  DISPERSION  VARIANCE,  S1GF 
125  IF ( X  .GT,  WTO)  IF ( T A  -  t F)  3 UO, 200, 200 
I f  (  T  A  .  GT.  TIE)  G)  to  !!■  0 


SIGTA  2 
SIGTA  3 
SIGTA  4 
SIGTA  5 
SIGTA  f 
SIGTA  7 
SIGTA  6 
SIGTA  9 
SIGTA  10 
SIGTA  11 
SIGTA  1? 
SIGTA  13 
SIGTA  14 
SIGTA  15 
SIGTA  16 
SIGTA  17 
SIGTA  18 
SIGTA  19 
SIGTA  20 
SIGTA  21 
SIGTA  22 
SIGTA  23 
SIGTA  24 
SIGTA  25 
SIGTA  26 
SIGTA  27 
SIGTA  28 
SIGTA  29 
SIGTA  30 
SIGTA  31 
SIGTA  32 
SIGTA  33 
SIGTA  34 
SIGTA  35 
SIGTA  36 
SIGTA  37 
SIGTA  36 
SIGTA  39 
SIGTA  40 
SIGTA  41 
SIGTA  42 
SIGTA  ^3 


o  o 


SIGE  =  ( SIGC1S  +  r  A*FPSCl)*-*3 
GO  TO  500 

150  StGE  =  3IGC2S  *TA> EPFC2 
GO  TO  500 

200  IF ( T A  .GT.  TLL)  GO  TO  25  0 

SIGE  =  (SIGC1L  +T1  *EPSC1) **1 
GO  TO  500 

250  SIGE  =  SIGC2L  +TA* £PSC? 

GO  TO  500 

300  SIGO  =  (SIGOP4-  (Tl  -  TO)  *SINT)*-*FX 

IFtTA.GT.  (1000.0  -  SI  GO ) / E  PSC 1)  GO  TO  35  0 
SIGE  =  (SIGO  +TA*rPSCi)‘*3 
GO  TO  500 

350  SIGE  =  3 . 0E6*SI GO  +TA*FPSG2-  2.0E9 
CORRECT  CLOUD  H0RI70NTA.  DISPERSION  VARIANCE  FOP  HIND  SHEAR  DISPERSION 
C  AND  RETURN  THE  STANDARD  DEVIATION 

5  00  SIGY  =  SORT  (S  IS  f  *■  <  SHEA  RC*  TA)  *•  2) 

RETURN 

END 


SIGTA  44 
SIGTA  45 
SIGTA  46 
SIGTA  47 
SIGTA  46 
SIGTA  49 
SIGTA  50 
SIGTA  51 
SIGTA  5? 
SIGTA  53 
SIGTA  54 
SIGTA  55 
SIGTA  56 
SIGTA  57 
SIGTA  56 
SIGTA  59 
SIGTA  60 
SIGTA  61 
SIGTA  62 


FUNCTION  ON AFl <  V ,  WIL'D,  SIF-E,  W,  I N  T  L  2 ) 

C 

COMPUTES  TOTAL  FALLOUT  FRACTION  (INTfGPATEO  IN  THE  CROSSWIND  OIFlCTION) 
DEPOSITED  AT  DISTANCE  Y  FRCH  C2  ALONG  THE  HOTLINE  (PER  6) 

H.  G.  NGRMiNT,  AT‘TO?P  H£  FI C  SCIENCE  ASSOCIA~ES  -  AUGUST  ISM 
C 

DATA  PI/3. 141592654/ 

I F  ( I  NT  L  2 )  10,10,  5  orr 

COMPUTE  INITIALIZATION  PAFAM.TEPS  ON  FIRST  PASS 
10  INTL2  =1 

WL  ~  A  LOG  ( W) 

IF  ( M  .GT.  10.0)  GO  TO  31) 

A  =  570.0 
B  =  0. 0176 
GO  TO  50 

30  A  -  -3179. d2  +  36  O'0 . 0  *HL 
B  =  0.  03 045*W**-  (-I.  GO 
50  IF ( W  .GT.  1.0E3)  GO  ^0  200 
IF  (W  .GT.  0.1)  GO  Tn  60 
ALPHA  “  1.06 
GO  TO  100 

60  ALPHA  =  +  'L  P 1 194  3 1*  WL 

100  IF  ( H  .GT.  1.0)  GU  "'0  150 
TM  -  3  0.  -*  W**<  0«-if-c6) 

GO  TO  TOO 

15U  TM  =  30.  *  H'MO.M  407) 

GO  to  .,0  0 
200  ALPHA  =  1.17 

IF  ( W  .  GT.  l.llt-O  GO  TO  300 
TM  =  5  3  3.61fe*W‘"(0-.l(273) 

GO  TO  -tOO 

300  TM  =  2497  1'3»WMT, 05116) 

400  Xl  =  TM*3QRT  (ALPHA/ (4.  0  -  ALPHA))  *  W^n 

COF  =  2, 9  *  i  IN (PI* ( 3, 0  -  A  L  PH A) /2, 0 ) / ( 3, 0  -  ALPHA)/PI**2 

PIX1  -  PI-Xl 

X1S  =  XJ.**2 

FX  IF  -  4.1*XIS**2 

fx is  --  -,.0‘xis 
txis  -  i.o-'xis 

VOX  =  WIN  )*W**(  0-  ’FFL'UJ/O.  66T-5 


DNAF1  2 
ON  AFl  3 
ON  AFl  4 
DNAF1  5 
DNAF1  6 
ON  AFl  7 
ON AFl  £ 
ON AF 1  9 

ON  AFl  in 
ON AFl  11 
DNAF1  12 
DNAF1  13 
DN  AFl  14 
ON AFl  15 
CN  AFl  If 
ON  AFl  17 
ON AFl  16 
ON  AFl  19 
ON  AFl  21) 
ON AFl  ?1 
ONAF1  2? 
DN  AFl  23 
ONAF1  24 
ON  AF 1  25 
ON AFl  26 
ON  AF 1  27 
0  N  AFl  26 
ON AFl  29 
9NAF1  30 
DN  AF  1  31 
ON AFl  3? 
ON  AFl  3.3 
ON AFl  34 
ONAFi  35 
ON  AFl  .36 
ONAFI  .37 
ONAFI  36 
ONAFI  3Q 
ONAFI  4P 
ONAFI  41 
ONAFI  42 
ONAFI  4.3 


"  4 


VTAli  =  2463)  *  HIND 

IF  (  W  -  3d,7866)  450',  AFO,  a?0 

U20  VATC  =  -1;  <i43*WIN9*F  XP(  lO.^VOF  +  WE*1  (  0  .  It  6176*  - 
1  fl.  660-4<44l-3*WL)  ) 

GO  TO  301 

450  l/ATC  =  -1.443*WIND*  14667.0*WJ‘*  (0.26203) 

Q  *  *  *  .  J *.»*»»«♦*»»*.*.**»**«•*■»■***•***»♦♦♦* 

COMPUTE  G ( X ) 

500  XS  = 

X2S  =  X 3  4-  S  IGF 
X3S  =  X2S  -  X1S 
X3F  =  X3S**2 
SRSIG  =  SORT(SIGE) 

XOSIG  =  X/3PSIG 
fXSXIS  =  XS*FX IS 

FXSX1F  =  XS*FX1F 
I F  (  X  3 F  :gt.  1v  0)  30  TO  i  10 
FI  =  X3S 
F2  -  X3F 
GO  TO  550 
510  FI  =  1.0/T5S 
F2  =  1,0 

FXSX1S  =  FXSXIS/XT 
FXSX1F  =  FX3X1F/X7F 

550  0NAF1  ~  SRSIG  *  COF*  <-X*  PIX1*(  <  3. 0*X2S  +  X1SJTT  +  FXEX1S)  - 

1  ( F  2  *  FXSXIS)*^?  4-  (X??*F2  -  rXSXin  *ALOG ( X2S/X 1 S )  + 

2  XOSIG*1  ( (  X2S  TX1D  *X2S*F1  4-  FX*X1F)MPI  +  ?.  1*  AT  AN  ( X  DS IG ) )  >  / 

3  (  (  F  2  +  3.  0*FXSX15>)  *  X3F  4  h,  O^XT^FX^XlF) 

CHECK  FOR  UPHINO  OR  DOWNWIND  CORRECTION 

IF (X) 600,000, 700 
CONFUTE  UPWItO  CORRECTION 

bOO  ONAF1  =  0NAF1  *  EX'Pt  P*X*(1,0  -  EXP(X/A))) 

RETURN 

COMPUTF  dornwino  correction 

700  UNAF1  =  OMAF1  ♦  EXT < - ( 6 L PG ( VOK^ON AF 1 )  4-  (  X/tfTAU)  “*  2)  * 

i  (i.o  -  -:xp(x/vato  ) ) 

e  0  U  RETURN 
END 


ON AF 1  44 
DNAF1  UF 
DN AF 1  AF 
DNAF1  A7 
ON  AF  1  AP 
ONAF1  AQ 
ONAF1  50 
DNAF1  51 
DN  AF 1  5? 
DNAF1  53 
ON AF 1  5A 
DNAFt  55 
ON  AF  1  66 
ONAF1  57 
DN  AF 1  56 
ON AF 1  59 
ONAF1  6 0 
ONAF1  61 
ON  AF 1  6? 
ONAF1  63 
CNAF1  6  A 
ONAF1  65 
ON AF 1  66 
PNAFl  67 
ON AF 1  66 
ONAF1  69 
DNAF1  70 
ONAF1  71 
ON  AF1  7? 
ONAF1  73 
DNAF1  7A 
DNAF1  75 
ONAF1  76 
DNAF1  77 
DNAF1  7B 
ON  AF  5  79 
DNAF1  60 
DNAF1  81 


» 


4 


- 


>• 


ft_ 


116 


'•  1 


O  O  O  O  C;  Cj  n  o  o 


1  YMIN,  YMAX,  OGX,  OC-V) 


H.  G.  N0R1.HT ,  ATMOSPHFF'.IC  SCIFNCF  ASSOCIATES 


1  R/  HR  CONTOUR  LEV7EL,  AND  S£TS  UP  the  MAR  ACCORDTNGL V, 
THE  HAP  IS  SPATIALLY  UNT IST0R"£0. 


DATA 
1  / 


NGRIOX 

IT 


c 


NGF.IDY 

44 


IH 

in 


IV 

6 


cw»*i  0.271P5  + 


TtNWNQ  -  10,0-WIN'? 

XON  =  (2.51E3  *  TINHND^* ( 0, 6834) ) 

1  A  LOG ( T£ N NND) )  ) 

IF (W  ,GT.  10.0)  GO  TO  100 

XUP  =  -3mO.O  *  (TrNWMP**  (0.1062)  )  *  <W**(0.15)) 

GO  TO  200 

100  XUP  =  -IPS. -5  *  (  T  JF  WHO*  *  ( 0  „  1 06  7)  )  M  ( H*+  (  0.  56  33)  ) 
COMPUTE  MAXIMUM  4n0  minIWM  HAP  C OOROI NA T FR 
?00  IF(SINA)  300,400,400) 

30  0  XMAX  =  XUP  *  SItIA 
XMIN  =  XCN  ♦  SIHA 
GC  TO  500 

L00  XMAX  =  XON  *  SIMA 
XMIN  ■=  XUP  *  SIMA 


600  IF (COS A) 600, TOO, 700 


600  VMAX  =  XU3  ♦ 
VMIN  =  XON  * 
GO  TO  300 
700  YMAX  =  XON 
Y M I  N  -  XUP 


COSA 

COSA 


COSA 
nnsA 

COMPUTE  GRID  INOREMEN’S  ANC  ADJUST  HAP  BOUNDARI FT 
LOO  XLNGH  —  XMAX  -  xmtti 
VLNGH  =  YHX  -  Vhiin 
IF (XLNGH  .GT,  Y  LNGH)  GO  TO 
P5  Y  -  YLNGH/NGP.20Y 
OGX  =  2, 0*OGY*IV/iH 
OG  =  ( 0, 567-VLNGH 
IF (OG) 1000, 1000, 8^0 
c  S  0  XHAX  =  XMAX  -*■  OG 
XMIN  =  XMIN  -  PG 
GO  TO  10 00 

,00  OGX  =  XLNGH/NGF.ID)’ 

OG  Y  =  OGX*IH/IV/?,iO 
OG  -  (0,667*XLNGH  -  YLNCH)/2.0 
IF (OG) 1000, 1000,9:0 
,50  YMAX  -  YH  AX  +  nc 


900 


XLNGH) /2. 0 


VMIM 


YMIN  -  PC- 


1000  XMIN  =  XMIN  -  OGX 
YMIN  =  YMIll  -  ICY 
RFTIIRN 
FN  U 


XMAX, 

SETMP 

? 

setmp 

3 

SFTMF 

4 

PATTERN 

sethp 

5 

SETMP 

F 

1981 

SETMP 

7 

SFTMP 

8 

CES  TO  THE 

SETMP 

9 

SETMP 

10 

SFTMP 

11 

SETMP 

12 

SETMP 

13 

SETUP 

14 

SETMP 

15 

LEVFL 
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DIMENSION  /MAP ( 20) ,  APSSA(tO),  H0LL(12),  OMAP(  5PQ0) 

DATA  INC/  13/,  XG/.,  YC  Z/P»  0,0.0/ 

format  c/itx, lores 

FORMAT  (1  X  ,  F 13.  n  ,  PX',  19 Ft  3> 

FORMAT (//15X,?5HTHE  CUf  U  T I TY  PRESENTED  IS) 
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026  RETURN  TRPL 
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044  MB  =  MA  T  RDL 

GO  TO  024  TRPL 
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SUBROUT  IN-:  ERROR  <PFOGPM,IP.RCR,ISOUT)  error  2 

;  T.  W.  SCHWENKE  ERROR  2 

Z  1  MARCH  1366  ERROR  *• 

:  ERROR  5 

;  *♦*♦**♦**■**< *-»-*+»*»-»*»l»*»»v*-**  •******•»**-*-»■("*“■»-*****■•****  **o ******•+*♦*£  RROR  6 

;  ERROR  7 

C  THIS  PROGRAM  WQITCS  A  GrNEPALT.IFn  FRP.OR  COMMENT  OF  THE  E  CL  LOWING  ERROR  6 

C  FORM  ON  TAPE  ISOUE  A  ED  TH£N  RETUP  NS  IE  THE  SIGN  OF  IRROR  IS  ERROR  9 

C  POSITIVE  OR  STOPS  IF  ITS  SIGN  IS  N-GA^IVS.  ERROR  10 

0  ERROR  11 

0  iRROR  SENSED  IN  F P.OC- HAM  <PROGRM)  At  OR  NEAR  STATEMENT  HUMBER  ERROR  1? 

C  (IRRCR),  PLEA'S!  REFER  TO  "HE  PROGRAM  LISTING,  ERROR  13 

C  ERROR  14 

C  PRIOR  TC  CALLING  -P.POR  THE  PARAMETEP.  PROGRH  MUST  BE  SET  ERROR  15 

Q  WITH  THr  RPO  NAME  OF  THE  CALLING  ERROR  IE 

C  PROGRAM  ANO  PApAMlTEP  IPROR  MUST  PE  SEf  WITH  THE  NUMBER  OF  THE  ERROR  17 

C  FORTRAN  STATEMENT  WHICH  PEST  IDENTIFIES  THF  FKRDP  CONDITION.  ERROR  IP 

C  ERROR  19 

C  ERROR  71 

1  FORMA'r(//?bH  ERP.O-'  SENSED  IN  cROGPAM  A6.30H  AT  OR  NEAR  STATEMENT  ERROR  7? 

1  NUMBER  IojAOH  .  PLEASE  REFER  TO  THE  PROGRAM  LISTING.)  ERROR  73 

C  ERROR  ?U 

ERROR  75 

Q,*******ERROp  26 

ERROR  27 

l'P.R-  I  A  IS  ( IRROR)  ERROR  76 

WRITE(ISCUT,l)Ppn'.'R(  ,  IPP.  ERROR  ? 9 

IF ( TOROS) 101, 100, 10  0  ERROP  30 

100  RETURN  ERROR  31 

101  STOP  ERROP.  37 

FN  n  ERROR  33 
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